{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "bd4214a3",
   "metadata": {},
   "source": [
    "# EoN库的使用：ER网络和BA网络上的流行病仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "87eef0f1",
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/usr/bin/python3\n",
    "# -*- coding: utf-8 -*-\n",
    "# Author ： 单哥的科研日常\n",
    "# 关注B站和公众号：单哥的科研日常，获取更多讲解教程"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bc0e50e5",
   "metadata": {},
   "source": [
    "### 实验环境：\n",
    "### Python版本==3.9.6, networkx==2.6.3, matplotlib==3.5.2, EoN==1.1 "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "68c6a16f",
   "metadata": {},
   "outputs": [],
   "source": [
    "import networkx as nx\n",
    "import matplotlib.pyplot as plt\n",
    "import EoN"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "01c13a3b",
   "metadata": {},
   "source": [
    "### ER均匀网络上的仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "fdb20433",
   "metadata": {},
   "outputs": [],
   "source": [
    "N = 10**4\n",
    "M = 4*N\n",
    "ER = nx.gnm_random_graph(N, M)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "13bfe6a2",
   "metadata": {},
   "source": [
    "#### SIS仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "3b2dcfc5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'individuals')"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEICAYAAAB1f3LfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAr3UlEQVR4nO3dd5xU1f3/8deHBSkiHSyAgoIFbMiCKEQRCyRRMUUl/lRsIRps0Wj0m1hiYmKNLdFIlNhQYokRE0UNQWwEWLAgoA8RFVZQuoBI2d3P749zJ7sLW2Z37szd8n4+HvdxZs5tn1mW+ey959xzzN0RERGprSZJByAiIvWbEomIiGREiURERDKiRCIiIhlRIhERkYw0TTqAXOvUqZP36NEj6TBEROqV2bNnr3T3zhWta3SJpEePHhQUFCQdhohIvWJmn1W2Tre2REQkI0okIiKSESUSERHJiBKJiIhkRIlEREQykrVEYmbjzWy5mb1fpq6Dmb1iZh9FZfsy6642s4Vm9qGZDS9T39/M5kbr7jYzi+qbm9nfovoZZtYjW59FREQql80rkoeAEdvUXQVMcffewJToPWbWBxgF9I32udfM8qJ97gPGAL2jJXXMc4E17t4LuAO4OWufREREKpW150jc/bUKrhJGAkOj1w8DrwK/iOonuvtm4BMzWwgMNLNPgTbuPh3AzB4BTgJejPa5PjrW08Afzcw8W+Piv/EGvPwy5OWVLk2alH9f2VLRdk2bQsuW0KpVadmmDbRtC82bZ+UjiIhkQ64fSNzZ3ZcBuPsyM+sS1XcF/ltmu8Kobmv0etv61D5LomMVmdlXQEdg5bYnNbMxhKsadt9999pFPn06/OY3tdu3ppo3DwmlXTvo0AE6doROncLrDh2gS5ewdO4MffpA+/bVHlJEJFvqypPtVkGdV1Ff1T7bV7qPA8YB5Ofn1+6K5YorwlJSAsXFpWU6S0XbFhXBpk2wcWPpsm4drF0LX30VlrVrYfVqWLoU3nsP1qyBDRu2j23vvaF3bxg+HI4/Hnr2rNVHFBGpjVwnki/NbNfoamRXYHlUXwh0L7NdN2BpVN+tgvqy+xSaWVOgLbA6m8ED4TZVkwQ7u23eDMuXw4oVoSwogKlTYeZM+Ne/4OKL4aijwnLooTBsWLiNJiKSJbn+RpwEjI5ejwaeK1M/KuqJ1ZPQqD4zug223swGRb21ztxmn9Sxfgj8J2vtI3VJ8+bQvTsccgiMGAG/+hVMmQJffgkffwzXXgvLlsF114UrlHbtwnZ33VXx1YyISIYsW9+9ZvYEoWG9E/AlcB3wD+BJYHdgMXCyu6+Otv8lcA5QBFzq7i9G9fmEHmAtCY3sF7m7m1kL4FGgH+FKZJS7L6ourvz8fG8UgzauXw+TJ8Orr4Zy0aLQpvLjH8M558BeeyUdoYjUI2Y2293zK1zXGP6IL6vRJJJtvfEG3HwzvPBCaLPp2xeefDI01ouIVKOqRKIn2xuLIUPg+edh8WK45BKYNw8OPBBuugka2R8TIhIvJZLGpmtXuPNO+OIL+MEP4Oqr4dxzlUxEpNaUSBqrnXeGiRPDcyl//StceGG45SUiUkNKJI2ZGaxcGZ6PuffekEx0ZSIiNaRE0tiZhUb4886D++4Lt7qUTKQC33wT/u6A8CtSUBB+fWqznHNOeP429atWXBx6r++2G/zkJzB3Ljz3HPzyl2H7AQPgrbfgs2iy16IiePZZuOGG0Cnxyithn33Kn2P//UNPeAidFvfYo/z67t1D2akTjBoV+p/88IcwZkw4r3tYNmyA2bPh0UfD9vfeG54ldoetW2HLFrjnHrj/fnjzzfAscUlJWL95MyxcCBMmhI6UEGK64ILa/dxatKh6fa9e4ef22GPbr7vzziz+crh7o1r69+/vUoHiYvdzzw3/d66/PuloGqTNm91XrHDfuNF94cLaH6e4OCwlJe6vv+7ep0/4Z1u0KBx73brwTwjuV1zh/sEHlR9ryRL3c85xf+IJ9/vvT311ammoywMP1P73Dihwr/h7tcLKhrwokVShuNj99NPDr8W4cUlHU6eVlIQf1+23hx9X8+bJf0kksVxzjfv48e4XXuj+9dfua9aEn03KmjXlf25ff+3+m9+kf/z773d//PHqt2vd2v2ss0I877zj/tRT7oMHl9+mZ0/3995zv+WWkDz79w/vn3oqJOETTnA/9lj3nXeu/Dw77RTKJk0y+7mNHOn+yivlf1bbWr/efeVK9+nT3SdNCkngmmvCf9FFi9y/+KL8/hs3lj/H+eeH9evXu599tnvfvu7z59f+d76qRKLnSKS8jRvDNf/q1TBtGhxxRNIR1RmffQY9eiQdxfYGDIBZs7av//Of4fzza3aszp3D0G4aVSc7li8P/Vvq489Xz5FI+lq1CjeUAUaObNTDqnz8cfhPn7rHXFkSOfzw8HhObf42Td1LnzAhfIGXXbd1axin85VXQhKr7BgzZ1Zc/5OflL5ON77ly+vnl1x90aVLw/z5KpHI9tq2hddfD6MPX3JJ0tHk1MKFIZemGi7XrCm/fvbs7b9833wzXMTVhkVjWJ92Guy6a/l1TZuGodKOOQZqO/tBSm3jE0mHEolUbMgQGDQIxo8v7frSgC1dGr7Ue/cOvZNSvve98PFTVw6HHJJcjCJ1lRKJVO7RR2GHHeCyy5KOJGs2bQoJpGvX0ro99ghNRMXF8Pe/wy67lF45iMj2lEikcr16hflNnn46/MnegJSUhOTQsmVp3dy54arj00/DpJNJTjsjUp/ov4pU7YILQnnLLcnGEaPPPoO8vNL3u+wSGrb33z+5mETqMyUSqdqee8L/+38wblzo0lPP/eMfpb2vBg4MVybLljXMnjQiuaJEItW7+urQAj1mTNKR1Jp7uJX1ve+F9/37w4wZavsQiYMSiVRvn33gsMPCID5r1yYdTY0tXFi+vWPNmjBOlIjEQ4lE0vPHP4by3nuTjaOGTjkldOkFGDw4DPbXrl2iIYk0OEokkp5DDoERI+COO8o/aFGHnXcePPVUeD1xYphtuGwju4jEQ4lE0veLX4RxxO+4I+lIqvWrX8GDD4bXK1fCqacmG49IQ6ZEIuk78kg4/nj4wx/CRAt11D/+ATfeGF5Pnw4dOyYajkiDp0Qi6TODsWNh1Sp4/vmko6nQqlWlPbMKCsIoLyKSXUokUjPHHhvGExk/PulIKrTffqG8777QxVdEsk+JRGomLw/OOgteegk+/zzpaMpZtgxWrAjTkdZ0Hg4RqT0lEqm50aPDI+FPPpl0JOXstlsoX3010TBEGh0lEqm53r2hXz94/PGkI/mf1atLXx96aHJxiDRGSiRSO6edFlqz581LOhIgTP4EpV1+RSR3lEikds48E5o1g4ceSjoSfvlLePvt0FvrnHOSjkak8TF3TzqGnMrPz/cCDbQUj/btw9hbqck9ErBlCzRvHl5v3Fh+fhERiY+ZzXb3/IrW6YpEau/440M5f35iIdx1Vyh/8hMlEZGkKJFI7f3ud6GcNCmR0xcXw5/+BEcdBX/+cyIhiAgJJRIz+5mZzTOz983sCTNrYWYdzOwVM/soKtuX2f5qM1toZh+a2fAy9f3NbG607m4zzS6RU927w0EHwc03J3L6iy8Osx2OHZvI6UUkkvNEYmZdgYuBfHffH8gDRgFXAVPcvTcwJXqPmfWJ1vcFRgD3mllqDNf7gDFA72gZkcOPIgCtWsFXX8GXX+b0tN98Uzqi/Ukn5fTUIrKNpG5tNQVamllToBWwFBgJPBytfxg4KXo9Epjo7pvd/RNgITDQzHYF2rj7dA89Bh4ps4/kynXXhfI//8npaVNJZOpUDQ0vkrScJxJ3/xy4DVgMLAO+cveXgZ3dfVm0zTKgS7RLV2BJmUMURnVdo9fb1m/HzMaYWYGZFaxYsSLOjyPHHBN6bOWwG/C6daF55rjjYOjQnJ1WRCqRxK2t9oSrjJ7AbsCOZnZ6VbtUUOdV1G9f6T7O3fPdPb9z5841DVmqkpcXJkR/+eVQ5sC114Yn2X/725ycTkSqkcStrWOAT9x9hbtvBf4OHA58Gd2uIiqXR9sXAt3L7N+NcCusMHq9bb3kWqqx/cMPs34qd5g8Gfr2hQEDsn46EUlDEolkMTDIzFpFvayOBhYAk4DR0Tajgeei15OAUWbW3Mx6EhrVZ0a3v9ab2aDoOGeW2Udy6bTTQvnMM1k/1TPPhHx12WVZP5WIpCmJNpIZwNPAHGBuFMM44CbgWDP7CDg2eo+7zwOeBOYDk4Gx7l4cHe4C4AFCA/zHwIu5+yTyP926hTndp0zJ+qluvz2UJ5+c9VOJSJo0RIrEY+xYePTRMGRKk+z8fbJ+PXTqBOedFx5EFJHc0RApkn2HHBK+6RcsyNopbrstjK31wx9m7RQiUgtKJBKPvn1D+eyzWTvFDTeEUl1+ReoWJRKJx6GHQpcu8MEHWTn8+++HcuDAxAYaFpFKKJFIPMxg2LDwhHsW2t1SzzvmoGOYiNSQEonE56ijYNky+OijWA9bVBR6a/XtGzqIiUjdokQi8TnyyFC+9lqsh506NZTnnRfrYUUkJkokEp+994add4Zp02I97JQpYVZfJRKRukmJROJjBkccERJJjO0kN98cEknr1rEdUkRipEQi8TriCFiyJMw4FYPp00PZoUMshxORLFAikXil2kliur318suhfOyxWA4nIlmgRCLx6ts3XD7ElEgKCmCvvUrzk4jUPUokEq8mTeBb34ql59bmzeGxlOHDY4hLRLJGiUTid+SR8PHH8PnnGR3mzTdh40YYMSKmuEQkK5RIJH7f+lYon3gio8NMnhx6ax11VAwxiUjWKJFI/A4+OJSvv57RYV56CYYMUbdfkbpOiUTi17QpjBoF//0vlJTU6hBLl8J77+m2lkh9oEQi2XHccbB8ea3nJ0l1+1VDu0jdp0Qi2XHYYaGcNatWu0+eDLvsAgceGGNMIpIVSiSSHXvvDW3awIwZNd61uBheeSVcjWjuEZG6T4lEsqNJkzALVS0Sybx5sHo1HH10FuISkdgpkUj2DBoUWsy//rpGu6WeZTziiCzEJCKxUyKR7Bk0KNynmj27Rru98Qa0awe7756dsEQkXkokkj2HHhrK1BC+afrb32DtWrWPiNQXSiSSPZ06hUb3t95Ke5dFi0I5bFiWYhKR2CmRSHYNHhwGzUrzwcRJk0L5l79kMSYRiZUSiWTXt74Fq1bB/PlpbX7PPWFIlD33zHJcIhIbJRLJriFDQvnii9VuWlISbm3tt1+WYxKRWCmRSHb16hXKK6+sdtOPPw7l2WdnMR4RiV3TpAOQBs4MDj8cvvii2k1TnbtSFzEiUj+kdUViZieb2U7R61+Z2d/N7JDshiYNxnHHwSefwPr1VW725pvQtm2YrVdE6o90b21d4+7rzWwIMBx4GLivtic1s3Zm9rSZfWBmC8zsMDPrYGavmNlHUdm+zPZXm9lCM/vQzIaXqe9vZnOjdXeb6cmDOmnAAHAPE7BX4c03w1iPTXTDVaReSfe/bHFUfhe4z92fA3bI4Lx3AZPdfV/gIGABcBUwxd17A1Oi95hZH2AU0BcYAdxrZnnRce4DxgC9o0WzV9RFAwaE8vbbK91kzZowxtbgwTmKSURik24i+dzM7gdOAV4ws+Y12LccM2sDHAE8CODuW9x9LTCScKVDVJ4UvR4JTHT3ze7+CbAQGGhmuwJt3H26uzvwSJl9pC7p3DmUVcxNkmofUSIRqX/STQanAC8BI6Iv/Q7AFbU8557ACuCvZva2mT1gZjsCO7v7MoCo7BJt3xVYUmb/wqiua/R62/rtmNkYMysws4IVK1bUMmzJSMeOoW+ve4Wr33gD8vLCgMEiUr9UmUiidosOQAvgVWBV9H4zUPUN78o1BQ4h3CLrB3xNdBursjAqqPMq6revdB/n7vnunt859dex5NZV0T/xkiUVrv7970OO2XHHHMYkIrGorvvvbEq/tLf98nbC1UVNFQKF7p6aqOJpQiL50sx2dfdl0W2r5WW2715m/27A0qi+WwX1UhelxoSfMWO7YX1To6eot5ZI/VTlFYm793T3PbcpU0utBrFw9y+AJWa2T1R1NDAfmASMjupGA89FrycBo8ysuZn1JDSqz4xuf603s0FRb60zy+wjdU2/ftCsGcyZs92qVNPJ5ZfnOCYRiUXaDyRG3XF7E25zAeDur9XyvBcBE8xsB2ARcDYhqT1pZucCi4GTo3PMM7MnCcmmCBjr7qleZBcADwEtgRejReqiZs1g61a4885wH6uMVEN7app3kYZq69atFBYWsmnTpqRDqVSLFi3o1q0bzZo1S3uftBKJmZ0HXEK4ffQOMAiYDtRqsG93fwfIr2BVhZOruvuNwI0V1BcA+9cmBklA585QQWeHqVPDsyO9eycQk0gOFRYWstNOO9GjRw/q4mNv7s6qVasoLCykZ8+eae+Xbq+tS4ABwGfufhTQj9DzSiR9l10WyrVry1U//nhoJ6mD/69EYrVp0yY6duxYJ5MIgJnRsWPHGl8xpZtINrn7puhEzd39A2CfavYRKe+gg0L57rv/q0pN5z56dAXbizRAdTWJpNQmvnQTSaGZtQP+AbxiZs+hHlJSU/nR3cwyQ6WkpnM/+eQE4hFphG688Ub69u3LgQceyMEHH8yMGTOq36kaabWRuPv3opfXm9lUoC0wOeOzS+PSuXN46nD8+P910UrllPyKWsxEJFbTp0/nn//8J3PmzKF58+asXLmSLVu2ZHzcdBvby3b8/yQqdyH0rhJJ3157lZst8frrQ7nzzsmEI9KYLFu2jE6dOtG8eXMAOnXqFMtx0+3++y9KH0hsAfQEPiQMpCiSvjPOgGuuga++grZtWb8+TK0r0uhceim88068xzz44NDFvhLHHXccN9xwA3vvvTfHHHMMp556KkceeWTGp02rjcTdD3D3A6OyNzAQeCPjs0vjc0g0jc277/6vJ/C11yYXjkhj0rp1a2bPns24cePo3Lkzp556Kg899FDGx63VDInuPsfMBmR8dml8+vUL5Tvv8N+vwrAphx6aYDwiSaniyiGb8vLyGDp0KEOHDuWAAw7g4Ycf5qyzzsromOm2kVxW5m0TwqCLeo5Eam6XXaBLF3j7baZ+CjvsUDpdiYhk14cffkiTJk3oHT39+84777DHHntkfNx0r0h2KvO6iNBm8kzGZ5fGxwz23x/mzeOOWaGRvWXLpIMSaRw2bNjARRddxNq1a2natCm9evVi3LhxGR833e6/v874TCIpvXuz+P4XAPj2txOORaQR6d+/P2+99Vbsx60ykZjZ81QyxweAu58Ye0TS8HXsyDRCT5HLLqtmWxGp86q7IrktKr9PeG7ksej9j4BPsxSTNHTnn88Vvwsji+63X8KxiEjGqkwk7j4NwMx+4+5HlFn1vJnVdgh5aey6dqWVfQYOTWvVb1BE6pJ0x9rqbGb/m8gqmmBKc9ZKrbg14asm7Rmz44SkQxGRGKT79+DPgFfNbFH0vgfwk6xEJA3e8uWwurgdfb+eASU/CpORiEi9lW6vrclm1hvYN6r6wN03Zy8sacgefzyUe/ExLFoEvXolG5CIZKTKPwXNbFhUfh/4LrBXtHw3qhOpsQ8+COVg3oQJur0lkkutszC4XXX3FFKjeZ1QwXJ87NFIo7BkCRx4oNOOr+Duu5MOR0QyVF2vreuil+e5e3EO4pFGYN48GDzY4KAz4OmnoahI3bdE6rF0//d+YmaTgb8B/3H3Sh9SFKnK6tWweDGMHQvs8V149FF4+20NuCWNTgKjyGdNut1l9gH+DYwlJJU/mtmQ7IUlDVXqP87BBwOHHx7exDDVp4gkJ91eW98ATwJPmll74C5gGpCXxdikAXr33VAefDDQuVt4c8MNcOGFSYUkkoiERpHPirQ78JvZkWZ2LzCHMEviKVmLShqsBx8MZZcuhJGAmzblfzNciUi9lFYiMbNPgEuB14H93f0Ud9cw8lJjn322zdS6v44Gll67NolwRCQG6Ta2H+Tu67IaiTR469bBhg1w441lKvPzQzl7Nhx9dCJxiTQmGzZsiP2Y1Q0jf6W73wL81sy2W+/uF8cekTRYs2eHMjVtO1CaSKZOVSIRqaequyJZEJWzsx2INHypHlvlEkmHDqG88Ub47W9zHZKIxKC6BxKfj8qHcxOONGQPPBDKLl22WTFkCLzxBriHBngRqVc0Q6LkzPz5law45ZSQSAoLoXv3nMYkkmvuTkVNBXVFbZ43r67X1m3A7cAnwDfAX6JlA/B+jc8mjVZJCbRoUcnUuocdFsrp03Mak0iutWjRglWrVtXqyzoX3J1Vq1bRokWLGu2X2AyJZpYHFACfu/vxZtaBMARLD8I0vqe4+5po26uBc4Fi4GJ3fymq7w88BLQEXgAu0fAtddO8ebBpExx0UAUrU5W33RauTkQaqG7dulFYWMiKOvzsVIsWLejWrVuN9km3+29nM9vT3RdBbDMkXkJozG8Tvb8KmOLuN5nZVdH7X5hZH2AU0BfYDfi3me0dDSJ5HzAG+C8hkYwAXswwLsmCWbNCOWhQBSubhfnbK7/3JdIwNGvWjJ49eyYdRuzSfbI9NUPiq2b2KjCV8IBirZhZN8L8Jg+UqR4JpBr1HwZOKlM/0d03u/snwEJgoJntCrRx9+nRVcgjZfaROuamm0JZ6RxWl14KxcWwZUuuQhKRmKSVSNx9MtCbcBVxCbBP6vZSLd0JXAmUlKnb2d2XRedbBqT69nQFlpTZrjCq6xq93rZ+O2Y2xswKzKygLl9SNmQffRRGQ6l0Vt3Bg8O9rzlzchqXiGSuJpNl9yfcXjoIONXMzqzNCc3seGC5u6f7bEpF3Ru8ivrtK93HuXu+u+d37pzpHTmpqXXRmAjf+lYVGw0eHMo338x6PCISr3TH2nqU0INrCDAgWvJrec7BwIlm9ikwERhmZo8BX0a3q4jK5dH2hUDZPqHdgKVRfbcK6qWOST3RftFFVWy0667Qs6cSiUg9lG5jez7QJ44eUe5+NXA1gJkNBX7u7qeb2a3AaOCmqHwu2mUS8LiZ/YHQ2N4bmOnuxWa23swGATOAM4F7Mo1P4pdqaK/yigTCVcnLL+vBRJF6Jt1bW+8Du2QzEEICOdbMPgKOjd7j7vMIc6HMByYDY8tM+3sBocF+IfAx6rFVJ82cGS42OnWqZsPBg2H5cvj445zEJSLxsHQuMsxsKnAwMBPYnKqvj0+25+fne0FBQdJhNCqpi4tqf9XmzoUDDwzT755+etbjEpH0mdlsd6+wSSPdW1vXxxeONCaFUb+6Sy9NY+N99w3luHFKJCL1SLpT7U7LdiDSML31VihPOy2NjVMPJi5YUPV2IlKnVNlGYmZvROV6M1tXZllvZproSqo1bVqYEbFfvzR3+O1vYdUq+OKLrMYlIvGpMpG4+5Co3Mnd25RZdnL3NlXtKwKhN++gQeFhxLSccEJoTPn977Mal4jEpyYPJIrUyPr1of089axhWg44IJR3352VmEQkfkokkjWzZoXh41OjxKfFDA4/HFq2zFpcIhIvJRLJmv/+N5QDB9Zwx5NOgm++Cc+UiEidp0QiWTN9eujR2759DXdMPQL/WkZT3ohIjiiRSFa4hyfaDz20Fjsfcgi0agWvvhp3WCKSBUokkhWLF4c7U7VKJDvsAMOGweTJscclIvFTIpGsSLWPDBhQywMMHx7G3Fq4MLaYRCQ7lEgkK/79b2jTppI52tMxYkQoX8pk/jQRyQUlEomdOzzwQOitlRr1pMZ69YI999TtLZF6QIlEYpcaKuvwwzM80IgRMHUqbN5c/bYikhglEond66+HMuMBfIcPh6+/1qyJInWcEonEbtq0MHNur14ZHuioo8K9MbWTiNRpSiQSK/fwHOERR8QwW+5OO8GQIWonEanjlEgkVosWweefw5FHxnTAAQPgvffg7bdjOqCIxE2JRGI1LZoCLbZEMmxYKM85J6YDikjclEgkVq+9Bp06wX77xXTA444L5UcfxXRAEYmbEonEatq0mNpHUszgwgtD48uWLTEdVETipEQisVm8GD79NMbbWilHHw0bN8LEiTEfWETioEQisUmN+n7EETEfOHV7609/ivnAIhIHJRKJzdSp0K5d6Wy5sWnVCs49N4xLv25dzAcXkUwpkUhs/vOf8AxhXl4WDn7WWaH86U+zcHARyYQSicRiwYLQPtKvX5ZOMGRIKCdMyNIJRKS2lEgkFjNmhHLo0CyeJDV41+zZWTyJiNSUEonEYu5caNECDjssiye5447QHfiZZ7J4EhGpKSUSicXMmXDwwdC0aRZP0qkT9OkDv/+9nikRqUOUSCRjq1bB9OlZvq2VcuGFoXz++RycTETSkfNEYmbdzWyqmS0ws3lmdklU38HMXjGzj6KyfZl9rjazhWb2oZkNL1Pf38zmRuvuNovteWqpgX/+E4qL4Qc/yMHJfvxj2G03uPvuHJxMRNKRxBVJEXC5u+8HDALGmlkf4Cpgirv3BqZE74nWjQL6AiOAe80s1cH0PmAM0DtaRuTyg0hw7rmhPOSQHJwsLw++/e3w9OPChTk4oYhUJ+eJxN2Xufuc6PV6YAHQFRgJPBxt9jBwUvR6JDDR3Te7+yfAQmCgme0KtHH36e7uwCNl9pEccQ9XIz17QpNc/TZdd11odL/++hydUESqkmgbiZn1APoBM4Cd3X0ZhGQDdIk26wosKbNbYVTXNXq9bX1F5xljZgVmVrBixYpYP0Nj9/77ofy//8vhSbt3hyuvhMcfD93FRCRRiSUSM2sNPANc6u5VjXtRUbuHV1G/faX7OHfPd/f8zp071zxYqdTPfx7KY47J8YmvvDJcDh14YChFJDGJJBIza0ZIIhPc/e9R9ZfR7SqicnlUXwh0L7N7N2BpVN+tgnrJoaXRT7xHjxyfuEOHMCowwKuv5vjkIlJWEr22DHgQWODufyizahIwOno9GniuTP0oM2tuZj0Jjeozo9tf681sUHTMM8vsIzmwciXMnw9XX51QAM8/Dx07wj33JBSAiEAyVySDgTOAYWb2TrR8B7gJONbMPgKOjd7j7vOAJ4H5wGRgrLsXR8e6AHiA0AD/MfBiTj9JI3fttVBSAiedlFAALVuGLmPPPgu/+11CQYiIeSO7v5yfn+8FBQVJh9EgdO8OhYVQVJSlEX/T8dVXYex6CN3HctZ1TKRxMbPZ7p5f0Tr9r5NaWb8e1q6FMWMSTCIAbdvCiSeG17/8ZYKBiDReSiRSK3/7G2zYAKNHV79t1j37LLRuDTfdpJGBRRKgRCK1MmEC7L13lkf7TVeTJvDuu+F1foVX3iKSRUokUmPLlsG0aXDaaeEB8zphzz1LX3/4YXJxiDRCSiRSYxMmhGcATz016Ui2kXrKPTX4l4jkhBKJ1Ig7jB8Phx8O++6bdDTb2H//UL75ZrhkEpGcUCKRGpkxI8zPfs45SUdSiQULQjl0aOhWJiJZp0QiNTJuHLRqBaecknQkldh33zAOF0D79rBpU7LxiDQCSiSStnffhb/+NTRB7LRT0tFU4eaboXfv8LplS03LK5JlSiSStltvDVcjv/510pGkoWzPrebNYdas5GIRaeCUSCQtn3wCEyfC+eeHO0Z1nlnoGdC2bXg/cCBcfHGyMYk0UEokkpaLLgpDWV12WdKR1NDatfDII+H1PfeEBKNbXSKxUiKRai1bBv/6V+jy27XCOSjruDPOKN/o3rw53HlnYuGINDRKJFKt228Po5Ck/rCvl5o3D7e6br01vP/Zz+DMM5ONSaSBUCKRKq1YAX/+M/zoR7DXXklHE4Of/xwefTS8fvTRcKvrN79JNiaRek6JRKr0l7/A1183sBHaTz89DF2ccu21IaGsXp1cTCL1mBKJVGr5crjttjA1+n77JR1NzHbcMdzqmjKltK5jx5BQJk5MLi6RekiJRCo1ZgysWdPAp0QfNiwklDfeKK370Y9CQjHT/CYiaVAikQp9+CE89xx85zsN8GqkIoMHh4Ry//3l6/PzS5OKWbhME5FylEhkO1u2lI7sO358srHk3JgxIaG4V/wI/847l08sl14KL74Yts9ESUnl64qKMju2SJYpkUg5xcWlsx5OnBi+Nxuta68tTSru4WGabd11V7hsa9KkfIKpbtljj/Lv8/Iq37ZZs8rX9ekDU6eGXhFl68eOhY0b4Zhjqo/lL3+BQYPC61tuyf3PWeo/d29US//+/V0qd8014VvziiuSjqSOW7TI/amn3Dt0KJtqGv7Svn3594sWuRcVuS9c6F5Skt7PbvVq97Fjtz/2gAHujzziPnJkad1777lv3er+059uv/3ZZ7v37eu+447uEya4v/OO+/r1IY7iYvdXXim//QsvhHVff+3+6afhXEVFYZ+iIvd//zsc47HH3A8+uHS/Zs3cDzssvH7xRfcTTnD/znfKH/uMM9z32sv98MMr/rkNH759LGvWpP/7tnWr+7p17p9/HuIvLt5+m5IS99Gj3Vu3Dj/DWbNK1y1d6n755e7LlqV/zm0ABV7J92oiX+ZJLkoklbvllvAbceaZ6X8nSBnffBP+gxcXhy+mxYtD/Zdfum/a5D55svuVV7o/8ID7rbe6b9gQftDb/ucuKan6H+DTT92nT3efMiV8gaW+nH73u/CFPmVK+S+tH/+49JglJe4rV4byyy/dzz8/bNOmTfaTkJbkl1tuqfWvtxJJmUWJpGI33RR+G77/fffNm5OORhKXSmRLl4bEWFDg/utfhyTl7j53rvuIEZl9qd16q/tHH4XjLV5cmtTAfdiwUG77F/4VV4TYZs1y79Mn1LVu7d6qVcXnuPzycJVxzjlVx9Kly/Z1c+a4T53qfv317oceWlq/227uzz7rPnNmKH/1K/cmTcK6kSPD1dHHH4ef2eTJ7ieeGPZ56KHwWV99Nd7kkO7Ss2dGvxJVJRIL6xuP/Px8LygoSDqMOuXBB0Mb83HHwfPPQ9OmSUck0sgUFZW2k33zTfjqb9Wq+v02bQrb52BIbjOb7e75Fa3TV0Yjd/zxoQ15+HB45hklEZFElP2P17Jl+vu1aBGWhOlro5EqKgpDTP3rX7D//jBpEuywQ9JRiUh9pO6/jcyWLWEQxk6d4IYb4NRTw8PbSiIiUlu6Imkkli6FO+6Ahx8OI/rm5YWrkBNOSDoyEanvlEgasPffD8+qvfQSTJ4c6kaOhPPOCw3reXnJxiciDYMSSQOyeTM8+yw8+STMmgWFhaF+jz1C8rjiigYyp4iI1Cn1PpGY2QjgLiAPeMDdb0o4pJzYuDGMH1hQEK463n8f5swJ02zsthsMHQoDB8JJJ8Huu4dehSIi2VCvE4mZ5QF/Ao4FCoFZZjbJ3ecnG9n2SkrCOFappagodAHfuDEs33wDX30VlrVrS1+vWlW6rF4NK1eGBLJxY+mxmzYNVx1nnhm68w4fHoZ+EhHJhXqdSICBwEJ3XwRgZhOBkUDsiWT8+DDddyoRbJsYqlqqGti1Ou3ahfmWOnaEzp1h773DQIpduoSeV336QP/+YUpyEZEk1PdE0hVYUuZ9IXDothuZ2RhgDMDuu+9eqxN16gQHHBAaqFNLkybl31e2VLZdq1ZhadkyLG3bli7t2sFOO6lBXETqvvqeSCq687/dmC/uPg4YB2GIlNqc6MQTwyIiIuXV9zvphUD3Mu+7AUsTikVEpFGq74lkFtDbzHqa2Q7AKGBSwjGJiDQq9frWlrsXmdmFwEuE7r/j3X1ewmGJiDQq9TqRALj7C8ALScchItJY1fdbWyIikjAlEhERyYgSiYiIZESJREREMtLo5mw3sxXAZ7XcvROwMsZw6gN95sZBn7lxyOQz7+HunSta0egSSSbMrMDd85OOI5f0mRsHfebGIVufWbe2REQkI0okIiKSESWSmhmXdAAJ0GduHPSZG4esfGa1kYiISEZ0RSIiIhlRIhERkYwokaTJzEaY2YdmttDMrko6nmwzs+5mNtXMFpjZPDO7JOmYcsHM8szsbTP7Z9Kx5IqZtTOzp83sg+jf+7CkY8omM/tZ9Dv9vpk9YWYtko4pbmY23syWm9n7Zeo6mNkrZvZRVLaP63xKJGkwszzgT8C3gT7Aj8ysT7JRZV0RcLm77wcMAsY2gs8McAmwIOkgcuwuYLK77wscRAP+/GbWFbgYyHf3/QnTT4xKNqqseAgYsU3dVcAUd+8NTInex0KJJD0DgYXuvsjdtwATgZEJx5RV7r7M3edEr9cTvly6JhtVdplZN+C7wANJx5IrZtYGOAJ4EMDdt7j72kSDyr6mQEszawq0ogHOqururwGrt6keCTwcvX4YOCmu8ymRpKcrsKTM+0Ia+JdqWWbWA+gHzEg4lGy7E7gSKEk4jlzaE1gB/DW6pfeAme2YdFDZ4u6fA7cBi4FlwFfu/nKyUeXMzu6+DMIfikCXuA6sRJIeq6CuUfSbNrPWwDPApe6+Lul4ssXMjgeWu/vspGPJsabAIcB97t4P+JoYb3nUNVG7wEigJ7AbsKOZnZ5sVPWfEkl6CoHuZd53owFeDm/LzJoRksgEd/970vFk2WDgRDP7lHDrcpiZPZZsSDlRCBS6e+pq82lCYmmojgE+cfcV7r4V+DtweMIx5cqXZrYrQFQuj+vASiTpmQX0NrOeZrYDoXFuUsIxZZWZGeG++QJ3/0PS8WSbu1/t7t3cvQfh3/c/7t7g/1J19y+AJWa2T1R1NDA/wZCybTEwyMxaRb/jR9OAOxdsYxIwOno9GngurgPX+znbc8Hdi8zsQuAlQi+P8e4+L+Gwsm0wcAYw18zeier+z91fSC4kyZKLgAnRH0mLgLMTjidr3H2GmT0NzCH0THybBjhUipk9AQwFOplZIXAdcBPwpJmdS0ioJ8d2Pg2RIiIimdCtLRERyYgSiYiIZESJREREMqJEIiIiGVEiERGRjCiRiIhIRpRIREQkI0okIgkysx+b2TvRUlLmdYMfTUAaDj2QKFIHRPNkvOXueyQdi0hN6YpEpG7YH5ibdBAitaFEIlI3HAC8X+1WInWQEolI3aArEqm3lEhE6gZdkUi9pcZ2kYSZWRNgLdDZ3TcnHI5IjemKRCR5vQizFCqJSL2kKxIREcmIrkhERCQjSiQiIpIRJRIREcmIEomIiGREiURERDKiRCIiIhlRIhERkYz8f9nAdTMLBycmAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "tau = 0.5          # transmission rate\n",
    "gamma = 0.2     # recovery rate per node\n",
    "t, S, I = EoN.fast_SIS(ER, tau=tau, gamma=gamma, tmax = 10)\n",
    "plt.plot(t, S, color = 'r', label='S')\n",
    "plt.plot(t, I, color = 'b', label='I')\n",
    "plt.legend(loc=0)\n",
    "plt.xlabel(\"$T$\")\n",
    "plt.ylabel(\"individuals\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "de302c28",
   "metadata": {},
   "source": [
    "#### SIR仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "9a15a8c6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'individuals')"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEICAYAAACuxNj9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA2nUlEQVR4nO3dd3xUVfr48c+TQhJCh9ASaRpQegkdpYkCKiBKsYKiYFnBXlb3t7u6KrvyXRUEXVQEhBVRkSpYWEBRigQQ6V0ILaHXAEnO748zMRGSTMlM7iR53q/XvO6dO+fe+4zCPJxzTxFjDEoppVReQpwOQCmlVPDTZKGUUsotTRZKKaXc0mShlFLKLU0WSiml3ApzOoBAqVSpkqlVq5bTYSilVKGSmJh42BgTc+nxIpssatWqxapVq5wOQymlChUR+S2n49oMpZRSyi1NFkoppdzSZKGUUsotTRZKKaXc0mShlFLKLU0WSiml3NJkoZRSyi1NFtllZMCHH8KMGU5HopRSQaXIDsrz2XvvQXIy9OwJkZFOR6OUUkFBaxbZhYTAyJGwZw9MnOh0NEopFTQ0WVyqSxdo0gT+8x/QVQSVUgrQZHE5EfjTn2DtWvj8c6ejUUqpoKDJIieDB0ODBjZpnD7tdDRKKeU4TRY5CQuzzVDJyTB5stPRKKWU47Q3VG7at7fPLiZOhEcecToapXKUYTI4n3ae0xdOczz1OKcvnOZc2jnOXTyX5/Z82nnSTTrGGDJMxh9ehsuP5XY8p/O9KZtZLvtnlx7Lfq1Lj+V1Xk4MOT+HNDk8n/SmbDBdG+Dws4eJDPNvb05NFnm56y549lnYtAmuucbpaFQRlWEyOHL2CIfOHOLw2cMcPnuYlDMppJxNIflMMilnU35/f+r8KVLTUklNS+Vc2jkupF/w6Z4RoRGEhoQSIiF/eAly2bEQCUHk8uO+lg0LCftDmeyfZR7PPJZTuUvL5PaZiOT43YVcjudQ3puywXTtUAnNsWx+SG6ZqbBLSEgw+V78KDkZatWC22/X5ijlFWMMh88eZsexHew5sYdDpw9x6Myh37fJZ5I5eu4ox1OPcyz1GGkZaTlep3xkeWKiY6gcXZmYkjGUiShDZFgkkWGRRIVF/b4fXSKacpHlKFWiFFFhUUSFR+W6jQiNyPUHSSkRSTTGJFx6XGsWealcGYYOhXHj4KWXoG5dpyNSQSQ9I50dx3aw7cg2dh/fza7ju9h5bCc7ju1g57GdnL7wx84RIRJCTMkYqpSqQuXoytQqV4vykeUpH1WeqqWqUiW6CjHRMVSMqkjl6MpUKlmJ8NBwh76dUn+kNQt3Dh6EK6+Efv10oF4xZIwh5WwKaw6sYeuRrRw4fYDfTvzGjqM72JCy4Q8JISI0gjrl63BlhSupU67O7/s1y9akaqmqVCxZkRDRPiUquGnNwldVq9qutBMmwFtvQblyDgekAiU1LZV1h9aRuD+R1QdWk3ggkY0pGzmffv73MqESSvXS1YmvGM+9je8loXoCV1e6mtrla1M5urImA1VkabLwxP3326aoadPgoYecjkb5wYX0CzYh7E8k8YB9bUjeQLpJB6BCVAWaV2vOn1r9iSvKXEGDyg1oWLmhJgRVbGkzlCeMsd1oo6JgxQr/XFMVqPSMdLYc2cLi3YtZtHsRX2//mlMXTgFQqWQlWlRrQYtqLWherTktqregZtma+hBYFUvaDJUfIrZ28cQTsGGDHd2tgt6h04f4cM2HfLfzO5YlLSM1LRWA2NKx9Kvfjx7xPWgd25q4MnGaGJRyQ5OFp+66C555xnah/ec/nY5G5WHtwbWM+mkU09ZPI8Nk0KRqE4a1GEaTKk3oWKsjtcvV1uSglJc0WXgqJsbOSDtjhp3GXH9sgsrZi2eZvmE6434ex8/7fyYyLJLhrYcztMVQrq50tdPhKVXoabLwRv/+8MADsHw5tG3rdDQK2HVsF2N/Hsv7q9/n5PmTxFeIZ0yPMQxoMICY6Binw1OqyNBk4Y3evW2y+O47TRYOuph+kblb5/KfxP/wzY5vEBFuu+Y2Hkp4iM61OmsTk1IBoMnCG5UqQYcO8PHHdkS3/igVqN+O/8b7q99nwpoJHDh9gLgycfy1418Z0nwIcWXinA5PqSJNk4W3Bg6061zs2AFXXeV0NEWeMYZZW2Yx+ZfJzNw8ExGhx1U9GNZiGD3iexAWon+ElSoI+jfNW9272+28eTBihLOxFHFLdi/h+YXPszxpOeUjy/Nc++d4uOXD1Chbw+nQlCp2NFl468oroX59mD1bk0WAHE89zhNfP8HEtROJKxPHB7d8wKCmg7QWoZSD9G+fL3r1glGj4PhxnSvKzz7f+DnD5g7jROoJXujwAn+57i9EhUc5HZZSxZ5OcuOLW26BtDRYsMDpSIqMY+eOMWjmIPp91o865evw84M/81rX1zRRKBUkNFn4onVr2zNqzhynIyn00jLSeOabZ4j9dyxT103lL9f9hZ/u/4lm1Zo5HZpSKhtthvJFaCjcfDPMnAnp6fa98trypOWMWDCClftWcm+Te3myzZM0qdrE6bCUUjnQmoWvunSxzyzWrXM6kkLp/376P9p92I5tR7Yx7bZpTOozSROFUkEsoMlCRJ4QkQ0isl5EPhGRSBGpICLfisg217Z8tvIviMh2EdkiIjdmO95CRH51fTZagmGIbpcudqtNUV5Jy0jjwdkP8vS3T9P3mr7sfnw3AxoOcDospZQbAUsWIhILDAcSjDENgVBgIPA8sNAYEw8sdL1HROq7Pm8AdAfGiUhm+867wFAg3vXqHqi4PRYbCwkJ8O23TkdSaBxPPU7nSZ35YM0HDG81nCl9p1AmoozTYSmlPBDoZqgwIEpEwoCSwH6gNzDJ9fkkoI9rvzcwzRhz3hizC9gOtBKRakAZY8wyY1dqmpztHGd16wbLlsGZM05HEvSSzyTTYUIHlu1dxke9P+LtHm8TGRbpdFhKKQ8FLFkYY/YBo4A9wAHghDHmG6CKMeaAq8wBoLLrlFhgb7ZLJLmOxbr2Lz1+GREZKiKrRGRVSkqKP79Oztq1sw+4V64M/L0KsS2Ht9D2w7bsPLaT+XfNZ3DTwU6HpJTyUiCbocpjawu1gepAtIjcndcpORwzeRy//KAx440xCcaYhJiYApieumNHCAvT8RZ5OHbuGDdMuYFT50+xePBiul3ZzemQlFI+CGQz1PXALmNMijHmIjADaAcccjUt4domu8onAVdkOz8O22yV5Nq/9LjzSpe2CWPmTLtOt/qDo+eO0nd6X/af2s+cO+bQKraV0yEppXwUyGSxB2gjIiVdvZe6ApuA2cAgV5lBwCzX/mxgoIhEiEht7IPsla6mqlMi0sZ1nXuzneO8m2+GrVth3z6nIwkqyWeS6TSxEz/u+ZEJvSbQOq610yEppfIhYIPyjDErRORzYDWQBqwBxgOlgOkiMgSbUPq5ym8QkenARlf5R40x6a7LPQxMBKKA+a5XcOjUyW6/+w4GD3YykqCx7cg2rv/4elLOpDBz4Ex6xvd0OiSlVD6JKaLNJwkJCWbVqlWBv5ExUL26TRqffBL4+wW51LRUWr3fiv2n9rPg7gUkVE9wOiSllBdEJNEYc9lfXJ3uI79EbBfa+fMhIwNCiu+g+NS0VPp/1p9fk3/lqzu/0kShVBFSfH/Z/KlzZzh8GDZvdjoSx5w6f4rOkzozd+tcxvYcS4/4Hk6HpJTyI00W/tCund3+9JOzcTgkPSOdvtP7snLfSqb0ncIjLR9xOiSllJ9psvCHunWhcmVYssTpSBzx8pKX+W7nd4zrOY47G93pdDhKqQDQZOEPItC+PSxf7nQkBe6NH9/g5e9f5t4m9zK0xVCnw1FKBYgmC39p1Qq2b4cjR5yOpMCsSFrB8wufp3e93rx303sEw2TASqnA0GThL23a2G0xqV0cPH2QvtP7Elcmjol9Juryp0oVcZos/KV1a4iKgq+/djqSgDPGcPeMuzl27hhfDviScpHlnA5JKRVgOs7CX6KioEMHWLTI6UgCbtIvk1i4ayFjeoyhebXmToejlCoAWrPwp44dYf36Iv3c4syFM7z0v5doXq25dpFVqhjRZOFPHTva7dKlzsYRQCOXjmTfqX2M6TGGENE/PkoVF/q33Z9atoSICPj+e6cjCYi1B9cy8seR3NXoLtpd0c7pcJRSBUiThT9FRNheUUVwcN6Rs0e46b83UTm6Mm/e+KbT4SilCpgmC39r1w7WroULF5yOxK+GzB7C4bOHmXfnPGKiC2AVQqVUUNFk4W+NGtl1uTdtcjoSv5m6biqztszilc6v0LRqU6fDUUo5QJOFv7VoYbcrVjgbh58cPnuY4QuG0yauDU+0ecLpcJRSDtFk4W/x8VClCvz4o9OR+MXby9/m2LljjL95POGh4U6Ho5RyiCYLfxOBtm3hhx+cjiTfUs6k8NaKt7it/m00qtLI6XCUUg7SZBEInTvDrl2we7fTkeTLXxf/ldMXTvPitS86HYpSymGaLAKhSxe7LcRdaL/b+R3vrnqXp9s+rQ+1lVKaLAKifn2IjoZVq5yOxCcnz59k0MxBxFeI5+XOLzsdjlIqCOhEgoEQEmJ7RRXSHlH/XPpP9p/az/Ihy3XqcaUUoDWLwGnb1g7OS011OhKv7Du5j7dXvE3/Bv1pHdfa6XCUUkFCk0WgtG4NFy/CmjVOR+KVZ759hrSMNF7v+rrToSilgogmi0Bp7fpX+cqVzsbhhWV7l/HJ+k94tv2z1Clfx+lwlFJBRJNFoFSvbgfnFZKaRYbJ4PGvH6d66eo82/5Zp8NRSgUZfcAdSE2awC+/OB2FR0avGM3KfSuZ0GsCpUqUcjocpVSQ0ZpFILVoYVfOC/KH3CfPn+TVH17l+jrXM7jpYKfDUUoFIU0WgdSiBaSlwa+/Oh1Jnkb9NIrDZw/zWpfXEBGnw1FKBSFNFoFUv77dbt7sbBx5OHn+JKNXjKbvNX1pGdvS6XCUUkFKk0UgxcdDVFRQP+R+b9V7nDh/ghc6vOB0KEqpIKbJIpDCwqBZs6DtPpualsq/l/2bbnW6kVA9welwlFJBTJNFoLVpY+eICsJlVj9a8xGHzhzSWoVSyq2AJgsRKScin4vIZhHZJCJtRaSCiHwrIttc2/LZyr8gIttFZIuI3JjteAsR+dX12WgpTE9h27WD8+eDrikqLSONN356g9axrelUq5PT4SilglygaxZvAwuMMVcDTYBNwPPAQmNMPLDQ9R4RqQ8MBBoA3YFxIhLqus67wFAg3vXqHuC4/adtW7tdtszZOC7x6fpP2XV8Fy90eEF7QCml3ApYshCRMsB1wIcAxpgLxpjjQG9gkqvYJKCPa783MM0Yc94YswvYDrQSkWpAGWPMMmOMASZnOyf4Va8ONWvCTz85HcnvMkwGI38cSf2Y+txS7xanw1FKFQKBrFnUAVKAj0RkjYh8ICLRQBVjzAEA17ayq3wssDfb+UmuY7Gu/UuPX0ZEhorIKhFZlZKS4t9vkx9t2wZVzWLe1nmsT17P8+2fJ0T0sZVSyr1A/lKEAc2Bd40xzYAzuJqccpFTW4jJ4/jlB40Zb4xJMMYkxMTEeBtv4LRrB0lJsHev+7IBZozh9aWvU7NsTQY2HOh0OEqpQsKjuaFEpB/22cMpEXkJmwT+YYxZncdpSUCSMSZzBaDPscnikIhUM8YccDUxJWcrf0W28+OA/a7jcTkcLzyyP7e44oq8ywbYot2LWJa0jHdvepfw0HBHY1GqKLp48SJJSUmkBvk0P5GRkcTFxREe7tnvgKcTCf7FGPOZiHQAbgRGYR8657o6jjHmoIjsFZF6xpgtQFdgo+s1CBjp2s5ynTIb+K+I/Buojn2QvdIYky4ip0SkDbACuBcY42HcwaFJEzs478cfoX9/R0MZs3IM0eHR3NP4HkfjUKqoSkpKonTp0tSqVStoO48YYzhy5AhJSUnUrl3bo3M8TRbpru1N2GalWSLyNw/OewyYKiIlgJ3Afdimr+kiMgTYA/RzBb9BRKZjk0ka8KgxJvO+DwMTgShgvutVeISH2/Utli51NIztR7czc/NM/tzhz0SXiHY0FqWKqtTU1KBOFAAiQsWKFfHm2a6nyWKfiPwHuB74p4hE4MHzDmPMWiCnocFdcyn/KvBqDsdXAQ09jDU4XXcd/OMfcOoUlC7tSAjvrHyHsJAwHm31qCP3V6q4COZEkcnbGD19wN0f+Bro7ur+WgF4xqs7FXetW0NGhh3N7YBT508xYc0E+jfoT/XS1R2JQSlVeOWZLFyjrSsAkcBi4Ijr/XnAmV+9wqpdOxCBJUscuf3kXyZz6sIpRrQe4cj9lVIF59VXX6VBgwY0btyYpk2bsmLFCvcnueGuGSqRrO6rl3ZjNdixFMoT5crZSQW//77Ab30h/QJvLn+TltVb0iq2VYHfXylVcJYtW8bcuXNZvXo1ERERHD58mAt+mJsuz2RhjPHsMbnyTKtWMHWqXRAprOBWtJ20dhI7ju1gdI/RBXZPpZQzDhw4QKVKlYiIiACgUqVKfrmux79Yrgn/4rFNUgAYYwr+n8mFWadO8N57sHq1TRwF5D+J/6Fxlcb0uKpHgd1TKQU8/jisXevfazZtCm+9levHN9xwAy+//DJ169bl+uuvZ8CAAXTs2DHft/XoAbeIPAB8j33I/XfX9m/5vntxk/k/7IcfCuyW6w6tI/FAIkOaDSkUPTSUUvlTqlQpEhMTGT9+PDExMQwYMICJEyfm+7qe1ixGAC2B5caYziJyNTZpKG9UrWpXz/v+e3jqqQK55YQ1EygRWoK7Gt1VIPdTSmWTRw0gkEJDQ+nUqROdOnWiUaNGTJo0icGDB+frmp52nU01xqQCiEiEMWYzUC9fdy6uOna0PaLS092XzafzaeeZsm4Kfa7uQ8WSFQN+P6WU87Zs2cK2bdt+f7927Vpq1qyZ7+t6WrNIEpFywEzgWxE5RmGbnylYXHstfPABrF9vpwEJoLlb53Lk3BHub3p/QO+jlAoep0+f5rHHHuP48eOEhYVx1VVXMX78+Hxf16NkYYy51bX7NxFZBJQFFuT77sVRV9fg9W+/DXiy+Hjdx1SOrsz1da4P6H2UUsGjRYsW/BSA9XM8fcBdI/MF7ALWAlX9Hk1xEBsLdesGfLxF0skk5mydw/1N7yc0JNT9CUoplQdPm6HmkTUoLxKoDWzBLoGqvNWhA8ycCcbYUd0B8H7i+xhjeLDFgwG5vlKqePGoZmGMaWSMaezaxgOtAGenUC3MWraEo0dhz56AXP7cxXOMXjmaW+rdQp3yOsheKZV/Pq2U51r0qKWfYyk+mje329V5rR3lu6+2fcXx1OM81uqxgFxfKVX8eLpS3pPZ3oZgV8oLokWuC5lGjSA01CaLW291X95LU3+dSuXoynSq1cnv11ZKFU+e1ixKZ3tFYJ9h9A5UUEVeVJRNGMuW+f3Se07sYdaWWQxuMpiwkIKbf0opVbR52nVWR2v7W7t28PHHdnBeqP96K436aRQAj7R8xG/XVEoVLqVKleL06dN+vWaeyUJE5mB7QeXIGNPLr9EUJ23awLhxsGkTNPTPIoCZCxzd3fhuapbL/4hNpZTK5K5mMcq17YsdVzHF9f4OYHeAYioeWre22+XL/ZYspq2fxpmLZ3g44WG/XE8ppTK5W89iCYCIvGKMuS7bR3NERKcnz4+rroLoaPj1V79d8v3V79MgpgGtY1v77ZpKKd85MEN5wHj6gDtGRH7vsC8itYGYwIRUTISE2BqFn/4kJe5P5Of9P/Ng8wd1KnKllN952l3mCWCxiOx0va8FDAtIRMVJmzYwfjxcvAjh4fm61ORfJhMRGsHgpoP9E5tSKt8cmqE8IDwdwb0Au0reCNernjHm60AGVixcey2cOwcrV+brMmkZaXy64VNuqnsTZSPL+ik4pZTKkmeyEJEurm1f4CbgStfrJtcxlR+dO9u5ob7OX979367/cejMIV3gSCkVMO6aoToC/wNuyeEzA8zwe0TFSYUKdp6oRYvydZmpv06lbERZesb39FNgSqnCzN9jLMB9b6i/unYfMMYEfmm34qhtW/vc4sIFKFHC69PPXTzHjE0zGNBgAJFhkQEIUCmlPO8NtUtExotIV9GuNv7Vvr19brFunU+nz906l9MXTnNnozv9HJhSSmXxNFnUA74DHsUmjndEpEPgwipGWrWy2xUrfDr9k/WfUCW6Ch1rdvRjUEop9Uee9oY6Z4yZbozpCzQDygBLAhpZcVGjBlSp4lOyOJ56nHnb5nFHwzt0NTylVEB5vJ6FiHQUkXHAauxqef0DFlVxImLHW/iwZu6MTTO4kH5Bm6CUUgHn6Rrcu4DHgR+AhsaY/saYLwIZWLFy7bWwYwccPOjVaVPWTeHK8leSUD0hQIEppZTlac2iiTHmVmPMJ8aYMwGNqDjq4Hr8s3ixx6esO7SORbsX8UDzB3R6D6XUH4SGhtK0aVMaNmzILbfcwvHjx/N9TXeD8p517f5DREZf+sr33ZWVkADlysHChR6fMmHNBEqElmBoi6GBi0spVShFRUWxdu1a1q9fT4UKFRg7dmy+r+muZrHJtU3M5eWWiISKyBoRmet6X0FEvhWRba5t+WxlXxCR7SKyRURuzHa8hYj86vpsdJHrvhsaapuilnjWZyAtI43//vpfetXrRYWoCgEOTilVmLVt25Z9+/bl+zruBuXNcW0n5eMeI7BJp4zr/fPAQmPMSBF53vX+ORGpDwwEGgDVge9EpK5rMOC7wFBgOfAV0B2Yn4+Ygk/HjjBnDuzfD9Wr51l0ye4lpJxNYWCDgQUUnFLKF48veJy1B9f69ZpNqzblre5veVQ2PT2dhQsXMmTIkHzf110z1BwRmZ3by93FRSQOO6fUB9kO9wYyk88koE+249OMMeeNMbuA7UArEakGlDHGLDPGGGBytnOKjo6ucRIe1C4+3/g50eHR9IjvEeCglFKF0blz52jatCkVK1bk6NGjdOvWLd/XDPRKeW8BzwKlsx2rYow5AGCMOSAilV3HY7E1h0xJrmMXXfuXHr+MiAzF1kCoUaOGB+EFkaZNoXRpmyzuuCPXYukZ6czcMpMbrryBkuElCy4+pZTXPK0B+FvmM4sTJ05w8803M3bsWIYPH56va+ZZszDGLHGtltfMGDPAGDPH9boTyHMEt4jcDCQbYzx6tgHk9BzC5HE8p3jHG2MSjDEJMTGFbG2msDDbK8pNzWLR7kUcPH2QgQ21CUoplbeyZcsyevRoRo0axcWLF/N1rUCulNce6CUiu4FpQBcRmQIccjUt4domu8onAVdkOz8O2O86HpfD8aKna1fYvBl++y3XIu+teo8KURW4ue7NBRiYUqqwatasGU2aNGHatGn5uo6nySJzpbzFIrIYWIQdpJcrY8wLxpg4Y0wt7IPr/xlj7gZmA4NcxQYBs1z7s4GBIhLhSkbxwEpXk9UpEWnj6gV1b7ZzipabXQlgfs7P7s9ePMv87fMZ2GCgNkEppXJ16RTlc+bM4Z577snXNT1aVtUYs0BE4oGrXYc2G2PO+3jPkcB0ERkC7AH6ue6xQUSmAxuBNODRbNOiPwxMBKKwvaCKVk+oTHXr2p5QixbBQw9d9vHsLbM5e/Es/Rr0cyA4pVRx5uka3AAtsGtvhwFNRARjzGRPTjTGLAYWu/aPAF1zKfcq8GoOx1cBDb2ItXASgRtugC+/hNRUiPzj+hTzt8+nUslKXFvjWocCVEoVV57ODfUxtmdUB6Cl66UTEgWA6defNSdqI1GRiNj8sXYtXEy/yPxt8+lau6vOMKtUkLO9/IObtzF6WrNIAOqbwvBfoBA7cgQq3dQD+OP4iWbNgDpL4N4UBjQY4EhsSinPREZGcuTIESpWrBi087YZYzhy5AiRkZ6vrulpsliPHWdxwJfAlHsPPAAffphHgfqfwflS9G3cHZO/HnBKqQCKi4sjKSmJlJQUp0PJU2RkJHFxce4LuniaLCoBG0VkJfD7g21jTC/vwlM5OXHij4nCTPsUBg6EpUuhfXuW/JBGpzkzYOvNkBZFzZqwe7dtolJKBZfw8HBq167tdBh+52my+FsggyjOjLETzgIMHgwffQSc7AHh4TBzJrRvz8XYxRB9mF5X9WM2sGcPhITYc5VSqiB4uqzqkpxegQ6uOHjxxaz9CRNcO2XK2AF6X34JxvDZhs+IDo9m2is9eOWVrPIi8OSTBRquUqqYcjeR4FLX9pSInMz2OiUiJwsmxKLr3Xfh9dftfnLyJc1KffrAjh2kr1/HrC2z6Bnfk6jwKF56CY4ezSr25pv2OkopFUju5obq4NqWNsaUyfYqbYwpk9e5yr1HHrHbkiXhsqmsetnHQctnvsOhM4e47Zrbfv+ofHnbBDVqVNZ1PvusAAJWShVbnk73ofxs5067rVsXzuS0UG21atCmDV9um0N4SHiO05E/9RQMcPWk7d8f5s4NXLxKqeJNk4VDrrzSbj/5JI9CffrwZflDdK3WnjIROVfkpk2DW2+1+7fcAuvW+TdOpZQCTRaO2LEja79589zL7eyWwM4K0ONYxTyvN2MGvz/4btIEXnjBD0EqpVQ2miwccNNNdjt6dN7lpp1ZAUCvb3KfsjzTSy/BiBF2f+TIrNqGUkr5gyaLArZuHWzZYvcfeyzvsp9v/Jw2xFHru1VZDzny8NZbsGaN3Z85E8aMyVeoSin1O00WBaxJE7t9/PG8y+04uoM1B9fQr+Vge+C///Xo+k2bwldf2f3hw2H58jyLK6WURzRZFKDx47P233wz77KfbbR9YW9v/yBce609OSPDo/v06AHPPGP327bVkd5KqfzTZFFADh2CYcPs/rZt7st/vvFzWse2pkbZGvDww7B3Lyxe7PH9/vUvKFvW7ofo/2WlVD7pz0gB6dLFbmvUgKuuyrvsgVMHSDyQSK96rnka+/SxI/G8HKp98GDW/t69Xp2qlFJ/oMmigJQqZbe/ue/YxLxt8wC4ua5rTe6oKDuH+YwZdhZBD0VGwsSJdr9GDS+CVUqpS2iyKADbtsHKlfDaa56Vn7Z+GjXL1qRR5UZZBx991G7/8x+v7j1oENSrZ/fdddVVSqncaLIoAHXr2u3Age7LHjt3jMW7F3N347v/uMpWzZp2iPb778OFC17df9kyux0xAr74wqtTlVIK0GQRcCdOZO17sh7K/O3zSTfpWU1Q2T3yCKSkwDvveBVD+fLw4492//bb4dQpr05XSilNFoHWvbvdvvWWZ+VnbJpBlegqtIptdfmH3brZdS5efx1SU72Ko107+OADu9+5s1enKqWUJotAyxwUN3y4+7InUk/w1bavuO2a2wiRHP7XiMCzz8Lhwz61Jw0ZYueiSkzU5iillHc0WQRQZtfV+HjP1sueuHYi59LOcV+z+3IvdP319iHIv/4F6elex/TNN3Z7++1/bCJTSqm8aLIIoA8/tNuxYz0rP23DNJpVbUZC9YTcC4WEwD/+YSeZynN+85xVrAjvvWf3n3jC69OVUsWUJosAeuklu80ckJeXpJNJLE9azu31b3df+PbbbX/YN97waS6PYcPs2t8ffQRz5nh9ulKqGNJkESCnT9ttxYoQGuq+/IxNMwA8SxYi8Pzztnbh4/J4/+//QfXqdvXW7OtrKKVUTjRZBMhzz9nt6697Vv6LTV/QqHIj6las69kJd95paxePP+51zyiAEiWypjC/8UavT1dKFTOaLAJk3Di7feAB92V3HdvFD7/9QL/6/Ty/QYkSdkj2zp1ZfWK91LevraDs2AG33ebTJZRSxYQmiwCYPt1ua9XyrBfUjE0zMBjuaXKPdzfq1g06dLBL450/73WcYJ+Vg512yst5CpVSxYgmCz8zBgYMsPuzZnl2zpytc2hUuRG1ytXy7mYi8Ne/wr59MGGCd+e6hIZmPbN45BHYvNmnyyilijhNFn723Xd2GxoKjRu7L7/v5D6+/+17br3ax0Wzu3a1w7Nfey3rqbqX6tSBVavs/jXXwLFjvoWilCq6ApYsROQKEVkkIptEZIOIjHAdryAi34rINte2fLZzXhCR7SKyRURuzHa8hYj86vpstIgnjTvOePllu123zrPyn6z/BIPh7sZ3+3ZDEdsMlZRkF/X2cVm8Fi2yppzq18+n8X5KqSIskDWLNOApY8w1QBvgURGpDzwPLDTGxAMLXe9xfTYQaAB0B8aJSGan03eBoUC869U9gHH77MIFWLrU7tev79k509ZPI6F6AvEV432/8bXXwp//bBev8LTtKwePPmpbsxYutF1rlVIqU8CShTHmgDFmtWv/FLAJiAV6A5NcxSYBfVz7vYFpxpjzxphdwHaglYhUA8oYY5YZYwwwOds5QSWzm+zjj3tWfmPKRhIPJHJnwzvzf/O//92ucPTXv+arWnDffXYOqdde83zkuVKq6CuQZxYiUgtoBqwAqhhjDoBNKEBlV7FYIPvin0muY7Gu/UuP53SfoSKySkRWpaSk+PU7eOJvf7Pbf/7Ts/Kfb/wcQRjY0IOFLtwJC7Mjutet83oK80uNGWOnNf/Tn2DmzPyHppQq/AKeLESkFPAF8Lgx5mReRXM4ZvI4fvlBY8YbYxKMMQkxMTHeB5sPR45k7Zco4b68MYYp66bQoUYHqpWu5p8g+vWDm26yIwI3bfL5MlFR8Msvdv/WW+HQIf+Ep5QqvAKaLEQkHJsophpjZrgOH3I1LeHaJruOJwFXZDs9DtjvOh6Xw/GgkrlkamZvKHcW717MtqPbeKC5B6P2PCViHzpERtoR3ufO+XypK66wCSMkBHr39vm5uVKqiAhkbygBPgQ2GWP+ne2j2cAg1/4gYFa24wNFJEJEamMfZK90NVWdEpE2rmvem+2coPFv1zf0ZNJAgHd+foeKURW9G7XticqVYfx4WLsWHnwQMjJ8vlTjxvDkk7BiBbRurQlDqeIskDWL9sA9QBcRWet69QRGAt1EZBvQzfUeY8wGYDqwEVgAPGqMyXxS+zDwAfah9w5gfgDj9lqyq27Utq1nI7aPpx5n1uZZDG46mKjwKP8H1L+/7cM7dar9tc+H116z3+vnnz1f7U8pVfSEBerCxpil5Py8AaBrLue8Cryaw/FVQEP/RedfU6bYrceTBm78gnSTTt9r+gYuqJdegqNH7S98zZo+L14RHg4//ABxcfDUU/bB9+DBfo1UKVUIBCxZFBfp6fD009CqFXTs6L58hsngjZ/eoEmVJrSNaxu4wERg1Cg7WO/JJ20b0hNPeFb1uUTmlCC9etmutWFhcLePYwiVUoWTTveRT/fdZ3+HB3rY+3XZ3mVsObKFx9s8TsAHooeGwn//a6eUfeopuP9+n59hlCwJ8+bZhPjgg7Bnj59jVUoFNU0W+fTpp3Y7bJhn5d9e8TalS5TmtmsKaE7w8HC7/OoLL9gR3sOG+ZwwIiLsUrEXLtiWLV3DW6niQ5NFPmzZYn84n3/e/svbnYOnDzJj0wyGtRhG6YjSgQ8wU3i4fVL94ot27YuePbOeynvpyiuzBh3WrGknvFVKFX2aLPLhgw9s+/2IEZ6Vn/zLZNJNun/HVnjjlVfs8OxFi+xDFk8HhVzi6aftuksnTthmqePH/RumUir4aLLw0ZEj9vlx165Qtar78sYYPlzzIR1qdKBepXqBDzAnInYOj6VLbZtSt252KT8ffu0fewwWLLAPvhs1gr173Z+jlCq8NFn4qK5rqex+Ho6pW7pnKVuPbGVIsyGBC8pTLVvaQXvPPWefY9Sv79NstTfeaJukkpLghhvg8GG/R6qUChKaLHxw4YIdwgC2g5EnxqwcQ+kSpf0/YttXUVF2HYwVK+yo7z597BJ/W7Z4dZlnn4XFi+0KezExsGZNQKJVSjlMk4UPMkcyDx3q2bCFeVvn8dnGzxjeejjRJaIDGpvXWrSww7P//ndbu2jY0GZAL9ZX7dgRvvrK7ick6EhvpYoiMUV0wp+EhASzKnOtUD9KTbX/KAc7IC/ETbo1xnDdxOvYc2IP2x7bRolQD6akdUpyMvzjH/bJfWqqnUHwnnugR4+sL52Ho0ftuJPZs+30VA8+WAAxK6X8SkQSjTEJlx7XmoWXMhc2evNN94kCYNHuRSzds5Rn2z0b3IkCbHPU6NHw22925b2FC+2Avrp1Yfp02/6WhwoV4IsvoFMnW+saONDtKUqpQkJrFl7KbHY6eRJKuxkqkZ6RTqsPWnHw9EF2DN9BZFik3+MJqIsXbffaJ5+0zVKlStlnG/362afbERE5npaaCp07w/LldsXXKVPsIn5KqeCnNQs/mDvXbkeOdJ8oAD7b+BmrD6zmjW5vFL5EAXYwX48edvW9WbPg9tvhyy9t81TNmrZ6lcOovMhIWLbMdrT6+Wdo104H7ylV2GnNwguZtYqjR+3sq3kxxtDy/ZacuXiG9Q+vJzQk1K+xOOb8eds89eabWYP66tSB9u3tQL8uXaBePTsvFfDjjzbflC1r803z5g7GrpRyK7eahc4666Hffsvad5coAJYnLSfxQCJje44tOokCbNNTz542AyxfbqsQ338P334LH3+cVaZhQ7jmGtrXqsX3z7Ti5ndupG3bcN57T7jvPme/glLKe1qz8FDHjvY38Zdf7Apy7nSZ1IVfk39lx/AdlIko47c4gpYxsGuXHXSxcSMkJtr3e/dCRgZHqMAAPmUh19O51M/MbPcvytSqYNdvrVrVLphRo4Z9lSrl9LdRqtjSmkU+/POfNlHcdptniWLV/lUs2r2If13/r+KRKMC20dWpY1/ZnT8P27dTcdMm5m9YyX2To5m6sy0VvvmEVeW60fT4+MuvVb68TRoVK9oZGitVss9IYmJsj62YGPsqV84+PCpd2qd1OpRSntOahRunTkEZ1+/9mTOezS5775f38uXmL9n35L7ikyy88N57tmfuuXPw8INpvHj/ASqe3Wvb+vbsyXodP27/ox88CAcO5H7BkBCbYCpUsEkkMhJKlLDNYSVL2jEiUVH2eOnS9gFKxYp2GxZmn69cui1RwpaPiLCvzOtlLxMaau+tiUoVIVqz8FGbNnY7ZoxniWLXsV18uuFThjQbookiFw89BLfcYpfYGD0ujDfHXMHw4Vfw4ovtqFw5l5MuXrSzN6ak2MGDKSm2//LJkzapHD1qXykpdnDHqVN2e/q07cubmmqz09mz/v9CmYnjwgVb28lMIp6+IiN/7xCglF9kThbqR1qzyMOPP0KHDnZmWU9m8zbG0OfTPizYvoDtj23nirJX5Ov+xcG6dXDrrbBzp31U0a0bjBvn2Uy+PklPt3OrJyfbWkt6OqSlXb49f/7y14UL9vPsZTP3z5yxCS001C4u5ekrPd0mMh8XpFIqR19+aWvDPtCahQ86dLDbmTM9K79g+wJmb5nNq11e1UThocaNYft223Hgz3+2f8a//NLW5B56yLb4+FVoqG2uqlDBzxdWqmjTQXm5yByAFx3tWeec5DPJDJ07lDrl6/B0u6cDG1wRIwJNm9rJCOfPt48dHnvMzjLi4/pMSik/02SRg6VLbZs6eLaojzGG+2fdT8qZFKbfPj3454AKYt27w6FDdooqY2yzVEKCXY1wxQp7TClV8DRZXCIjA554wu7v3evZALxxP49j3rZ5jLphFC2qtwhsgMWAiK1ZbN4MY8fa59qjR9vOBqVK2f8/iYmaOJQqSJosLjFlCqxaBZMm2XFi7mw7so2nv32anvE9ebTlo4EPsBiJiIBHHrFj+5KT7XiXatVsAklIgGuugVdftSv1KaUCS3tDZZOWBldeaXviLFvmfgryi+kXuW7idWw5vIX1j6yneunq+YhYeeroUTsV+tSpsGSJPda+vR068a9/2amplFK+0VlnPRAWBl9/bdf+8WRRo4fmPsTypOWM7TlWE0UBqlDBLqy0eDFs3QqvvAL799tFl66+2s5nOHmyNlMp5U+aLC5x9dXQqJH7cmNWjmHC2gk82+5Z7mh0R+ADUzmKj4eXXrLjNPbvh//7PzvubtAgqF0bhg2zPaoOH3Y6UqUKN22G8sHi3Yvp9nE3esb3ZEb/GUVrVtkiICPDToA7eTIsWpRVw4iLs+M6WrWy/yho21YXZVLqUrk1Q2my8NKWw1to82EbqpWqxrIhyygbWdbv91D+c/IkrFwJa9bYgX+rVsGWLVmfx8ZCs2Z2nEfTpnZlv1ynHFGqGNAR3H6wIXkD3ad2JzwknHl3ztNEUQiUKQPXX29fmU6ftgljyRKbRNassYMB09Pt53FxtimyZk27X6GC7fhQp459H1kIFz1UKr80WXjo6+1fc8cXdxAiIXxzzzfULl/b6ZCUj0qVghYt7CtTaiqsXm17wa1ZA+vX20GAR49efn7mDOnx8XDVVXabuR8bq3MCqqJJk4UbxhjeXvE2T3z9BHUr1mX2wNnUq6R9M4uayEi7Vni7dn88nppqBwVu3mwHae7da8d1HDoE27bBggV2jsFM0dH2uUhsbNbyG9mX4IiJsbWTMmV0ZnNVuBSaZCEi3YG3gVDgA2PMyEDfc+HOhbz4vxdZsW8Fver14pPbPqFkuAfzlKsiIzLS/vDHxub8eUaGTR7bttkJETdutM9GNmywAwmPHMn5vJAQW8MpW9Y2c5UpY5NKnTpQpYod61O1qt0vV86+oqM1wSjnFIoH3CISCmwFugFJwM/AHcaYjbmdk58H3In7E/nz//7MNzu+IbZ0LM+1f45HWz1KiGhPY+WdtDTbbTc52W5TUuy6TidO2CU3TpywCWXHDtvl9+DBP9ZUsgsNtUmldOms9ZwufZUsmbXeU0772deCymlpDRHPj+e3bOYxkayXcl5hf8DdCthujNkJICLTgN5ArsnCF8YYek3rxdytcyldojQju45keOvhRIVH+fM2qhgJC8uqJXjCGNuD6+BB+0pOtms7HT9uE8uJE/bzc+ey1nM6d85+nrm2U+b27NnCtUxGZsK4NJHklFw8+Sy35OPNcX9cw4lrr1nj97WPCk2yiAWyz/+aBLS+tJCIDAWGAtTwoQO9iBBfIZ6/d/o7I1qP0N5OqsCJ2KapsmXzP22JMXa9puzJI3tCyciwZS5djymnY/4um55u32cey/7Zpcd8/Sy3/yaeHvembDBdGwJTSyssySKnr37ZfyZjzHhgPNhmKF9u9O8b/+3LaUoFHZGsJcTLlXM6GlXYFZZG+CQg+9JzccB+h2JRSqlip7Aki5+BeBGpLSIlgIHAbIdjUkqpYqNQNEMZY9JE5E/A19iusxOMMRscDksppYqNQpEsAIwxXwFfOR2HUkoVR4WlGUoppZSDNFkopZRyS5OFUkoptzRZKKWUcqtQzA3lCxFJAX7z8fRKQFFYiFO/R3DR7xFcisr3AP9+l5rGmJhLDxbZZJEfIrIqp4m0Chv9HsFFv0dwKSrfAwrmu2gzlFJKKbc0WSillHJLk0XOxjsdgJ/o9wgu+j2CS1H5HlAA30WfWSillHJLaxZKKaXc0mShlFLKLU0W2YhIdxHZIiLbReR5p+PxlYhcISKLRGSTiGwQkRFOx+QrEQkVkTUiMtfpWPJDRMqJyOcistn1/6Wt0zH5QkSecP2ZWi8in4hIpNMxeUJEJohIsoisz3asgoh8KyLbXNvyTsboiVy+xxuuP1frRORLESkXiHtrsnARkVBgLNADqA/cISL1nY3KZ2nAU8aYa4A2wKOF+LuMADY5HYQfvA0sMMZcDTShEH4nEYkFhgMJxpiG2OUCBjoblccmAt0vOfY8sNAYEw8sdL0PdhO5/Ht8CzQ0xjQGtgIvBOLGmiyytAK2G2N2GmMuANOA3g7H5BNjzAFjzGrX/insD1Oss1F5T0TigJuAD5yOJT9EpAxwHfAhgDHmgjHmuKNB+S4MiBKRMKAkhWTFSmPM98DRSw73Bia59icBfQoyJl/k9D2MMd8YY9Jcb5djVxL1O00WWWKBvdneJ1EIf2AvJSK1gGbACodD8cVbwLNAhsNx5FcdIAX4yNWk9oGIRDsdlLeMMfuAUcAe4ABwwhjzjbNR5UsVY8wBsP/AAio7HI8/3A/MD8SFNVlkkRyOFep+xSJSCvgCeNwYc9LpeLwhIjcDycaYRKdj8YMwoDnwrjGmGXCGwtHk8QeuNv3eQG2gOhAtInc7G5XKJCIvYpugpwbi+possiQBV2R7H0chqWLnRETCsYliqjFmhtPx+KA90EtEdmObBLuIyBRnQ/JZEpBkjMms3X2OTR6FzfXALmNMijHmIjADaOdwTPlxSESqAbi2yQ7H4zMRGQTcDNxlAjR4TpNFlp+BeBGpLSIlsA/uZjsck09ERLDt45uMMf92Oh5fGGNeMMbEGWNqYf9f/M8YUyj/FWuMOQjsFZF6rkNdgY0OhuSrPUAbESnp+jPWlUL4oD6b2cAg1/4gYJaDsfhMRLoDzwG9jDFnA3UfTRYurgdEfwK+xv4FmG6M2eBsVD5rD9yD/df4Wterp9NBFXOPAVNFZB3QFHjN2XC856oZfQ6sBn7F/n4UiikzROQTYBlQT0SSRGQIMBLoJiLbgG6u90Etl+/xDlAa+Nb1d/29gNxbp/tQSinljtYslFJKuaXJQimllFuaLJRSSrmlyUIppZRbmiyUUkq5pclCKaWUW5oslFJKuaXJQqkCICIPZhsgmZFtv1COsFfFjw7KU6oAudaE+MkYU9PpWJTyhtYslCpYDbFTZShVqGiyUKpgNQLWuy2lVJDRZKFUwdKahSqUNFkoVbC0ZqEKJX3ArVQBEZEQ4DgQY4w573A4SnlFaxZKFZyrsCvmaaJQhY7WLJRSSrmlNQullFJuabJQSinlliYLpZRSbmmyUEop5ZYmC6WUUm5pslBKKeWWJgullFJu/X8Zr8+bVuI+IgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "tau = 0.5           # transmission rate\n",
    "gamma = 1.0    # recovery rate\n",
    "rho = 0.05      # random fraction initially infected\n",
    "\n",
    "t, S, I, R = EoN.fast_SIR(ER, tau, gamma, rho=rho, tmax = 20)\n",
    "plt.plot(t, S, color = 'r', label='S')\n",
    "plt.plot(t, I, color = 'b', label='I')\n",
    "plt.plot(t, R, color = 'g', label='R')\n",
    "plt.legend(loc=0)\n",
    "plt.xlabel(\"$T$\")\n",
    "plt.ylabel(\"individuals\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7b531928",
   "metadata": {},
   "source": [
    "### BA异质网络上的仿真"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9780edd0",
   "metadata": {},
   "source": [
    "#### SIS仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "ce3cd2fb",
   "metadata": {},
   "outputs": [],
   "source": [
    "N=10**4\n",
    "BA=nx.barabasi_albert_graph(N, 4) #create a barabasi-albert graph"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "17f1266f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'individuals')"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEICAYAAAB1f3LfAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAAul0lEQVR4nO3deZwV1Z338c8PGmi2Zt+kwUZBCDuyiKIRAxETjWgyKs+Tic5oYh7HaHTMopNJYhZfmknixCzmGZJoSDRRHzXR4BKViIaIKKtsIptKQytNA9LI2vR5/vhVyxUauPRdqm/19/163VfdPreWX3Xfrl9VnVPnWAgBERGRhmoWdwAiIlLYlEhERCQjSiQiIpIRJRIREcmIEomIiGSkKO4A8q1r166hrKws7jBERArKggULtoQQutX3WZNLJGVlZcyfPz/uMERECoqZvXWkz3RrS0REMqJEIiIiGVEiERGRjCiRiIhIRpRIREQkIzlLJGZ2j5ltNrNlKWWdzexZM1sdTTulfHaLma0xs1VmNiWlfLSZLY0++6mZWVTeyswejMrnmVlZrvZFRESOLJdXJL8Fzjuk7GZgVghhADAr+hkzGwxMA4ZEy9xtZs2jZX4JXA0MiF5167wK2BZC6A/8N/CDnO2JiIgcUc6eIwkhvFjPVcJUYGL0fgYwG/h6VP5ACGEvsN7M1gDjzOxNoCSEMBfAzH4HXAQ8FS1za7Suh4Gfm5mFXPWLP2cOPPMMFBX5q0WLg+/TeRUXQ+fO0K2bv1q0yEmYIiL5lu8HEnuEECoAQggVZtY9Ku8NvJwyX3lUtj96f2h53TIbonXVmNl7QBdgy6EbNbOr8asa+vbt27DI586F732vYcvWp2tXOOUUGDv24Kt/f2imaisRKSyN5cl2q6csHKX8aMscXhjCdGA6wJgxYxp2xfLVr8JXvgK1tVBTk/5r/36f7t4NW7dCZSW88w5UVMDy5TB9Otx1l29j4EC45hq46ipo165BYYqI5Fu+E8m7ZtYruhrpBWyOysuBPinzlQKbovLSespTlyk3syKgA7A1l8FjBs2b+6tVq+yss6YGVqzwK57p0+GGG+C73/X3n/lMdrYhsQjBX2b+Zz6eu5kh+DI7dkCnTodfqO7bB3v3+nnNpk0wfz6cfLJ/NTt1gupqaNsWdu3ysgMHoGNHv8PaoQNUVUGflP84s4PrXb0aFi+G9euhZ09o3drnX7IE+vXzdWzbBt27w6BBvr3OnX3eAwfgjTd8fbt2eQxmfv7Up4/P07693+0Nwc+vyst9W0VF/lnHjnDSSb6uVq18+QMH4P33Pb62bT2mQ+3f7/taXX3w91Vc7LG/+ab/fOCA3wzYvBnefdfLunXz9716QcuW/ntv1Qp27vTyDh08rlatYO1aGDDAY9iwAcrK/HdSW+t/39S/d/PmPn37bf9brVnj83bs6Nts2dLn37wZXn/df0c1Nb7+9u2hTRuP/6WXfB9OOsnPNdu08XU0JvlOJI8DVwB3RNPHUsr/YGZ3AifgleqvhBAOmFm1mY0H5gGXAz87ZF1zgX8C/paz+pFcKiqC4cP99cUvel3MTTfBP/0T3HorfOtbB//LJavqvi1bt8J778GWLX5Q69DBDxSdO/s/fk2NH6y3boUuXfyC8vXXoaTEl6uo8J+3bIFVq6B3bz/o7tzp6y8u9gNJ3fYmTPBtVFYe3HZpKWzc6AfD7ds/HGeLFn5Q2rfPD1jZ0rKlr7Numi9mB38Xx1J3vrZ374fLS0pgzx7/O4Vw+O8snw7dn/btPZnlUu/e/p3p3t2/a23bevlbb8HVV/vva948/7xlS0/8f/87PPAAXHZZ9uPJWSIxsz/iFetdzawc+DaeQB4ys6uAt4FLAEIIy83sIWAFUANcG0I4EK3qGrwFWGu8kv2pqPw3wO+jivmteKuvwnfmmTB7tn8bbr3Vj1p3361kchxC8LuGHTvCs8/CunWeADZuhKVL/cxw7Vo/q929O/vbP3AALr7Yz4LPPtvPpMG3XVnpZ7K7d/uZ5fDhflBcscKTVqdOfmbcoYN/vmmT3+VcudLPSk87zafduvnZcIsWvm+dO/vBY9cu/6qsWuXJoUMHvxKoqvJ5Kyp8vzt18vhat/bl3n0XZszwhPXNb3r8ZWVe/tZbvr3hwz3u6uqDiXDlSt+/HTt83TU1vq89e3pi7d7d96FZM9/+tm2eYDt29OWbN/crhLPP9jhmz/bk0ayZn2Nt3er71KWLv+/d2/+ezZr5Pp1wgn8+ezZMm+bLVlVBjx6waJEfbE86yWMbMeLgdmtqPLZ9+/yKoHNnj6WoyL8fO3Z4Qhg71uOurvb4Nmzw5Xbv9mXfe8+3t3at71fPnl42YIDva5s2MHSoz9+p08Fltm/3RFhU5NPJk329Bw4c/H5UVPiJTQhw6qm+3927++crV/qVX3W1X4m+847/TcaNg1/8wk9IwPe/Uye/ygSPKResEE/iMzFmzJhQEL3/hgBf/zr88Ifw05/CddfFHVGjtn8/3HcffOc7fuCrT8uWfitm0CA/0Lz1lpdNmOAHg379/KAyfrxfYWzb5snnhBP8n7G21g+spaUwZIgfsHv18oNby5b53V+RI6lL5uDJMVvnoGa2IIQwpr7PGktluxzKDO64w49oN9zgpyQTJsQdVaNRWwsvvuhtIObP9zPKrVENWb9+Pv3kJ/2Af/HF/s/Vs6efAaZj4sSchC2Sc+l+x7O6zfxvUtLWrBn87nd+TX7JJQdvzDcxVVV+m2LDBrjySi9Lva8/YYJXQp5/Plx0kVpQi+SbEklj17EjPPyw3/z87GfhL3+JO6K8+etf4fe/h/vv/3D5sGFwzjlewXj99X6lISLxUSIpBGPHek3iAw94LerAgXFHlDP79nnr58cf98rpOjffDJMmwemnH2yhIiKNg24CFIq77vJ2pD/8YdyR5MwLL8CoUXDbbd5a5fbbvXI7BH8/ebKSiEhjpERSKLp39wqC3//e21UmRAh+9dGtm1dw79rlzVA3bPCrkO7dj7kKEYmZEkkhueYav/fz+9/HHUlW/OUv3sZ+6lTPjRdeCMuWweWX68pDpJAokRSSoUP9ibSf/9yfXCpg99zjiWPFCrjzTu9G4rHHlEBECpESSaG57jp/xHXu3LgjabCHHvJ+KcF348YbP9zvk4gUFiWSQvOpT3m/FvfdF3ckDfLDH3oDtLPO8u4dxo+POyIRyZQSSaEpKfFKhT/9Kbs9+OXYgQPeC//Xvgaf/jQ8+aR6yhdJCiWSQjRlincKtWxZ3JGkpbbWE8iPfwzXXgsPPqgkIpIkeiCxEE2a5NPnnvMuWRuxffs83Dlz4Atf8HYCIpIsuiIpRH36+NPtzz0XdyRHVVvrbQPmzIEf/AD+53/ijkhEckGJpFBNmuTd3+ZzRKLjdPvtPtDjTTf5rS0NqSKSTEokhWryZB9RaN68uCOp1zPP+ABJ06YlulcXEUGJpHBNnOj9pTfC21uvvebtAU480W9n6UpEJNmUSApVp04wZkyjSyR19SIATz3VJIdPEWlylEgK2eTJfmtrx464I/nA9OledXPvvT6krYgknxJJIZs40Z/0e/nluCMBYOVKr1ifOBGuuCLuaEQkX5RICtlpp/kAzbNmxR0J1dUweLAPmXL//aoXEWlKlEgKWUkJnHFGo0gk3/62T++9F044Id5YRCS/lEgK3VlnwaJFsHt3bCG88463zrrsMu8aXkSaFiWSQjdypDeVirHfrR/8wJ+L/P73YwtBRGKkRFLoRo70aUyJpLoafvMbvxrp3z+WEEQkZkokha6szGu4Fy+OZfM/+5knk+uvj2XzItIIKJEUuqIiOPNMmD0775veutVvZ118MYwbl/fNi0gjoUSSBGedBUuXwvbted3sNdd4Hf+3vpXXzYpII6NEkgQTJkAIee3AsbbWx14vKjpYTSMiTZMSSRKMHetPAOYxkTz5pE8LdOh4EckiJZIkKCmBU06B+fPzsrkQ4LbbvJ7/05/OyyZFpBGLJZGY2Y1mttzMlpnZH82s2Mw6m9mzZrY6mnZKmf8WM1tjZqvMbEpK+WgzWxp99lOzJtwxx/jxMHeuH+VzbOZM797rq1+FFi1yvjkRaeTynkjMrDdwPTAmhDAUaA5MA24GZoUQBgCzop8xs8HR50OA84C7zax5tLpfAlcDA6LXeXnclcZl3DjYsgXKy3O+qS9+0af/+q8535SIFIC4bm0VAa3NrAhoA2wCpgIzos9nABdF76cCD4QQ9oYQ1gNrgHFm1gsoCSHMDSEE4HcpyzQ9dTXeOX6e5L33oLISLr0UWrfO6aZEpEDkPZGEEDYCPwLeBiqA90IIzwA9QggV0TwVQPdokd7AhpRVlEdlvaP3h5YfxsyuNrP5Zja/srIym7vTeAwb5hXuOU4kDz4INTXwla/kdDMiUkDiuLXVCb/K6AecALQ1s38+2iL1lIWjlB9eGML0EMKYEMKYbt26HW/IhaF9e++jJMeJ5De/gaFDfXBGERGI59bWZGB9CKEyhLAfeBQ4A3g3ul1FNN0czV8O9ElZvhS/FVYevT+0vOkaOTKniWTZMnjlFbjqKo03IiIHxZFI3gbGm1mbqJXVJGAl8DhQN67eFcBj0fvHgWlm1srM+uGV6q9Et7+qzWx8tJ7LU5ZpmkaOhHXrvCIjB+6911tpffazOVm9iBSoonxvMIQwz8weBhYCNcAiYDrQDnjIzK7Ck80l0fzLzewhYEU0/7UhhAPR6q4Bfgu0Bp6KXk1XXYX7a695tylZtG8f3HmnjzeS1LuDItIweU8kACGEbwPfPqR4L351Ut/8twG31VM+Hxia9QALVWrLrSwnkj//2aef+ERWVysiCaAn25OkVy+/XFiyJOurfvBB6NTJ60dERFIpkSSJWU4q3HfuhEcfhdNO05PsInI4JZKkGTHCm1ft35+1Vf7tbz796EeztkoRSRAlkqQZORL27oVVq7K2yiee8MdUbropa6sUkQRRIkmaLHeVEoInknPPhZYts7JKEUkYJZKkGTgQWrXKWiJZsgQ2boTzz8/K6kQkgZRIkqaoyPvdylIimTnTp5/8ZFZWJyIJpESSRHUtt7IwNsnMmd5DfY8eGa9KRBJKiSSJRo6EqirYlFnXY5s3e99auq0lIkejRJJEWapwf+opv6i54IKMIxKRBFMiSaLhw32ahUTSsyeMGpV5SCKSXEokSdS+PZx8ckaJpKYGnnkGpkxRl/EicnRKJEmVYVcp//gHbNum21oicmxKJEk1ciSsWQPV1Q1a/Otf9+mUKdkLSUSSSYkkqVLHJmmAV1/1afv22QlHRJJLiSSpMmi5tWMH1NbCd76T1YhEJKGUSJKqd2/o0qVBieSJJ3xaVpbViEQkoZRIkqpubJIGDHK1aJFPP/OZ7IYkIsmkRJJkI0fC0qXelvc4zJrlY4+0bZubsEQkWZRIkmzkSNizB954I+1F3n0XFi6EyZNzF5aIJIsSSZI1oML9rrt8qv61RCRdSiRJNnCgj0Z1HImkbta6HCQicixKJEnWogUMHXpcFe4rV8JFF0EzfTNEJE06XCTdwIFp15Fs2ABvvgnnnJPbkEQkWZRIkm7wYM8OO3cec9Y//MGno0fnNiQRSRYlkqQrLfXpCy8cc9bVq306ZkwO4xGRxFEiSbrTTvPp2rXHnHX1ahg7Flq1ynFMIpIoSiRJN3AglJTAqlVHnW33bnj5ZTj77DzFJSKJoUSSdM2awYgR/pThUcydC/v2qaJdRI5fWonEzC4xs/bR+/80s0fN7NTchiZZM26cd6C1f/8RZ/n73717rjPPzGNcIpII6V6RfDOEUG1mZwJTgBnALxu6UTPraGYPm9nrZrbSzE43s85m9qyZrY6mnVLmv8XM1pjZKjObklI+2syWRp/91EyDwtZr1CjYu/eozYAXLTp4F0xE5Hikm0gORNPzgV+GEB4DWmaw3buAp0MIg4ARwErgZmBWCGEAMCv6GTMbDEwDhgDnAXebWfNoPb8ErgYGRK/zMogpuQYN8unrr9f7cQg+kNWoUXmMSUQSI91EstHM/ge4FHjSzFodx7IfYmYlwEeB3wCEEPaFELYDU/ErHaLpRdH7qcADIYS9IYT1wBpgnJn1AkpCCHNDCAH4XcoykuojH/HpvffW+/GGDbBpE5xxRh5jEpHESDcZXAr8FTgvOuh3Br7awG2eBFQC95rZIjP7tZm1BXqEECoAomn3aP7ewIaU5cujst7R+0PLD2NmV5vZfDObX1lZ2cCwC1ibNj5dv77ej+fO9enpp+cpHhFJlKMmkqjeojNQDMwGqqKf9wLzG7jNIuBU/BbZKOB9ottYRwqjnrJwlPLDC0OYHkIYE0IY061bt+ONNxkuucSbZdVj7lxo3RqGD89zTCKSCEXH+HwBBw/ahx68A351cbzKgfIQwrzo54fxRPKumfUKIVREt602p8zfJ2X5UmBTVF5aT7nUZ9AgeOQRf2CkdesPffTSS/4gYosWMcUmIgXtqFckIYR+IYSTDpnWvRqSRAghvANsMLOBUdEkYAXwOHBFVHYF8Fj0/nFgmpm1MrN+eKX6K9Htr2ozGx+11ro8ZRk51LBhUFt7WIX77t3edbxua4lIQx3riuQDUXPcAfhtLgBCCC82cLvXAfebWUtgHfCveFJ7yMyuAt4GLom2sdzMHsKTTQ1wbQihrhXZNcBvgdbAU9FL6jN0qE+XLftQ86yFC/3xkvHjY4pLpAnZv38/5eXl7NmzJ+5Qjqi4uJjS0lJaHMctirQSiZl9HvgyfvtoMTAemAt87PjDhBDCYqC+rgEnHWH+24Db6imfDwxtSAxNTv/+PsjVkiXwuc99UDxnjk/1IKJI7pWXl9O+fXvKyspojI+9hRCoqqqivLycfv36pb1cuq22vgyMBd4KIZwDjMJbXkmhaNGi3q5SXnsN+vSBrl1jikukCdmzZw9dunRplEkEwMzo0qXLcV8xpZtI9oQQ9kQbahVCeB0YeIxlpLEZNsxvbaVYsABOVWc3InnTWJNInYbEl24iKTezjsCfgWfN7DHUQqrwDBsGlZXw7rsAVFd7rykayEqk6bjtttsYMmQIw4cPZ+TIkcybN+/YCx1DWnUkIYSLo7e3mtnzQAfg6Yy3Lvk1bJhPly2DHj1YssS7R9EViUjTMHfuXGbOnMnChQtp1aoVW7ZsYd8Rni87HulWtvdN+bHu8eieeOsqKRR1LbeWLoVJkz6oLlEiEWkaKioq6Nq1K62i0eu6ZqlyNN3mv09w8IHEYqAfsArvSFEKRY8e0K2bJxK8fqRnT+jVK+a4RJqiG27wh7iyaeRI+MlPjvjxueeey3e/+11OOeUUJk+ezGWXXcbZWRjNLq06khDCsBDC8Gg6ABgHzMl465J/w4d7E2C8AZfqR0Sajnbt2rFgwQKmT59Ot27duOyyy/jtb3+b8XrTfiAxVQhhoZmNzXjrkn8jRsAvfsGuHTWsWFHExRcfexERyYGjXDnkUvPmzZk4cSITJ05k2LBhzJgxg3/5l3/JaJ3p1pH8e8qPzfBOF/UcSSGKBrlaOvMtamtPVv2ISBOyatUqmjVrxoABAwBYvHgxJ554YsbrTfeKpH3K+xq8zuSRjLcu+Re13FowazuginaRpmTnzp1cd911bN++naKiIvr378/06dMzXm+6zX+/k/GWpHEYOBCaN2fhgkCXLv5Uu4g0DaNHj+all17K+nqPmkjM7C8cYYwPgBDChVmPSHKruBiGDGHh+o6MHg+N/CFbESkAx2q19SPgx/izI7uBX0WvncCyoywnjdjeEeNYtqMvp4464jmCiEjajnpFEkJ4AcDMvhdC+GjKR38xs4Z2IS8xW9ZjEvtpyan9qoAucYcjIgUu3b62upnZBwNZRQNMNdExawvfwmbeg/+ptjjeQEQkEdJttXUjMNvM1kU/lwFfzElEknMLq06kA9s5adMcjjAEjIhI2tJttfW0mQ0ABkVFr4cQ9uYuLMmlhUtbcGqb17FFC489s4jIMRyr1dbHQgh/M7NPH/LRyWZGCOHRHMYmOVBT4z2kfOmkKu9sS0SalHbt2rFz586srvNYVyRnA38DPlXPZwFQIikwa9bA3r0wfGQz+ONGH5ukR4+4wxKRAnasVlvfjt5+PoRwIA/xSI7VDZA4ZGJ3+COwaBGcd16sMYlIYUu3sn29mT0NPAj8LYSgBxAK1KJFUFQEQ6b29+YSCxcqkYjEIIZe5HMm3ea/A4HngGvxpPJzMzszd2FJrixZAoMGQXGPDtC/Px+MbiUi0kDpttraDTwEPGRmnYC7gBeA5jmMTXJgyRL4YBybU0+FV16JNR6RpiqmXuRzIt0rEszsbDO7G1iIj5J4ac6ikpzYvh3Kyw+OuMvo0fDmm7B1a4xRiUihSyuRmNl64Abg78DQEMKlIQR1I19g6irao57kD/Yhv2hRLPGISDKkW9k+IoSwI6eRSM5FQ7UfTCSjRvl04UKYpCfcRZqCbD9DAsd+IPFrIYT/Ar5v9fQ3HkK4PusRSc4sXQolJSljkHTpAieeqAcTRSQjx7oiWRlNdaRJgGXLvH7kQ+cEo0crkYhIRo71QOJfoumM/IQjuVJb61Uhl19+yAcjRsCf/gQ7d0K7drHEJiKFTSMkNhFvvum5oq5a5AOnnw4hwD/+AVOmxBGaSJMSQqC+qoLGoiHPm2uExCairsXWB01/65x2mt/revnlvMck0tQUFxdTVVXVoIN1PoQQqKqqori4+LiWi22ERDNrDswHNoYQLjCzzngXLGXAm8ClIYRt0by3AFcBB4DrQwh/jcpHA78FWgNPAl9W9y31q2uxNXjwIR+UlMCQITB3bt5jEmlqSktLKS8vp7KyMu5Qjqi4uJjS0tLjWibd5r/dzOykEMI6yNoIiV/GK/NLop9vBmaFEO4ws5ujn79uZoOBacAQ4ATgOTM7JepE8pfA1cDLeCI5D3gqw7gS6fXXoW9fzxuHOeMMePBBr0hplvYzqiJynFq0aEG/fv3iDiPr0j1q1I2QONvMZgPP4w8oNoiZlQLnA79OKZ4K1FXqzwAuSil/IISwN4SwHlgDjDOzXkBJCGFudBXyu5Rl5BD33QctWx7hw9NPh/fe82wjInKc4hoh8SfA14D2KWU9QggV0fYqzKx7VN4bv+KoUx6V7Y/eH1p+GDO7Gr9yoW/fvhmEXZj27fPpiSceYYbTT/fp3Ln13PsSETm647mPMRq/vTQCuMzMDm1ImhYzuwDYHEJI9+GF+po3hKOUH14YwvQQwpgQwphu3TK9I1d43njDp1deeYQZTjkFOneGl17KW0wikhxpXZGY2e+Bk4HFeIU3+EH7dw3Y5gTgQjP7JN75Y4mZ3Qe8a2a9oquRXsDmaP5yoE/K8qXApqi8tJ5yOUTdHauPfOQIM5jB2LF6MFFEGiTdK5IxwIQQwr+FEK6LXg3qHiWEcEsIoTSEUIZXov8thPDPwOPAFdFsVwCPRe8fB6aZWauokn8A8Ep0G6zazMabN8q+PGUZSbFkCTRv7uOQHNHIkbBiBezena+wRCQh0k0ky4CeuQwEuAP4uJmtBj4e/UwIYTk+FsoK4Gng2pRhf6/BK+zXAGtRi616vfKKPz/SuvVRZjrzTNi/H+bNy1tcIpIM6Tb/7QqsMLNXgA8q2TN9sj2EMBuYHb2vAurtgjaEcBtwWz3l84FDH7GTFCF41yif+tQxZjzzTG/6O3s2TJyYh8hEJCnSTSS35jIIyZ2NG6GyMqXr+CPp2NH7T5k9Ow9RiUiSpNv894VcByK5MWuWT3umc2Py1FPhV7+CAwe8UkVEJA1HrSMxsznRtNrMdqS8qs1MA10VgM1R27dzz01j5o4dfbpkSa7CEZEEOmoiCSGcGU3bhxBKUl7tQwj1dbYhjczy5dCrlz8mckyf/7xPFy7MaUwikizqWCnhVq2CgQPTnHnAAOjeHZ5/PqcxiUiyKJEkWG2t9/p7zIr2OmbeRcof/uDNvURE0qBEkmBr18L77/uzhmmryzqrV+ciJBFJICWSBFu82KfHlUhuvNGnTzyR5WhEJKmUSBJs8WIoKjrODn379fNKlWefzVVYIpIwSiQJtmSJd9R4nKNmwqRJ8MILB/ufFxE5CiWSBFu8+Dhva9UZOxZ27YInn8xyRCKSREokCVVZ6d2jNCiRfOITPl2xIpshiUhCKZEkVN3D6Q1KJD16+PTWW7MUjYgkmRJJQtW12BoxooErGDLEu5Xfm8mIyiLSFCiRJNTixVBaCl26NHAFd9zhU/UGLCLHoESSUA2uaK8zeTK0aQMzZ2YpIhFJKiWSBNq928dpzyiRFBfDOefAUxp0UkSOTokkgZYv9yFFMkok4K231q6FN97IRlgiklBKJAnUoK5R6nP++T7V7S0ROQolkgRasgTat/feTjJSVuatt9TvlogchRJJAi1e7M1+m2Xjr3vBBfDii/Dee1lYmYgkkRJJwtTW+hVJg58fOdQFF0BNjTpxFJEjUiJJmPXrobo6C/UjdcaPh06dVE8iIkekRJIwWator1NU5K23nnzSL3dERA6hRJIwixdD8+ZeR54155/vvUDef38WVyoiSaFEkjCLF8OgQdC6dRZXOmWKT++7L4srFZGkUCJJmIy7RqlPly5w9tmwZo1ub4nIYZRIEqSqCsrLc5BIAC6/HNatg8cfz8HKRaSQKZEkSEZjkBzLZz7jU9WTiMghlEgSJOMxSI6mQwe48kp4+mnYuTMHGxCRQpX3RGJmfczseTNbaWbLzezLUXlnM3vWzFZH004py9xiZmvMbJWZTUkpH21mS6PPfmpmlu/9aUwWL4YTToBu3XK0gSuv9CTy//5fjjYgIoUojiuSGuCmEMJHgPHAtWY2GLgZmBVCGADMin4m+mwaMAQ4D7jbzJpH6/olcDUwIHqdl88daWxyUtGe6owzYOBAuOeeHG5ERApN3hNJCKEihLAwel8NrAR6A1OBGdFsM4CLovdTgQdCCHtDCOuBNcA4M+sFlIQQ5oYQAvC7lGWanD17YOXKHCcSM78qmTNHXcuLyAdirSMxszJgFDAP6BFCqABPNkD3aLbewIaUxcqjst7R+0PL69vO1WY238zmV1ZWZnUfGosVK7xLrJwmEoDPfc6feLz33hxvSEQKRWyJxMzaAY8AN4QQdhxt1nrKwlHKDy8MYXoIYUwIYUy3nFUgxCvrXaMcSa9e/lzJHXd45hKRJi+WRGJmLfAkcn8I4dGo+N3odhXRdHNUXg70SVm8FNgUlZfWU94kvfqqj0Fy8sl52Ngtt/j0kUfysDERaeziaLVlwG+AlSGEO1M+ehy4Inp/BfBYSvk0M2tlZv3wSvVXottf1WY2Plrn5SnLNDlz53pHvVkZg+RYrrsOTjkFbr8dQr0XgSLShMRxRTIB+BzwMTNbHL0+CdwBfNzMVgMfj34mhLAceAhYATwNXBtCOBCt6xrg13gF/FrgqbzuSSOxZ4+P0z56dJ422Lw5/Md/+BOQzzyTp42KSGNloYmdUY4ZMybMnz8/7jCy6qWXYMIE+NOf4KKL8rTRvXuhuBi6dvWegUUk0cxsQQhhTH2f6cn2BHjpJZ+OH5/HjbZqBR/7GGzZ4n1wiUiTpUSSAHPmQP/+0LNnnjc8Y4ZXyvzoR3nesIg0JkokBa62Fp5/HiZOjGHjpaVw7bXwf/+vHlAUacKUSArcxo2wYweMqffOZR584xs+HO/Pfx5TACISNyWSAlfXdfzgwTEF0KMHTJvm/W9t3nzs+UUkcZRIClxdIslJ1/Hp+sY3vA3y978fYxAiEhclkgL3j3/ARz4CJSUxBjFwIFx2GfzsZzBvXoyBiEgclEgK2J498MILcM45cUfCwZZbn/0sVFfHG4uI5JUSSQGbPx927YJzz407Erwzx9mzYf16+Ld/izsaEckjJZIC9uKLPj3jjHjj+MDZZ3t9yX33qesUkSZEiaSAPfecV7I3qp7xb7nFK22mTIGtW+OORkTyQImkQO3a5RXtkyfHHckhWreGX/zC3597LuzbF288IpJzSiQFas4cP0Y3ukQCXvv/s5/BggXwrW/FHY2I5FhR3AFIwzz3HLRoAWedFXckR/ClL8Frr8EPfgBlZfB//k/cEYlIjiiRFKjnnvNK9rZt447kKH76U3j7bU8qvXrB1KlxRyQiOaBbWwVoyxYfo33SpLgjOYbiYnjwQejb1wdK+fGP445IRHJAiaQAPf+8j3DbKOtHDtWhAyxcCBdeCF/5ysGKeBFJDCWSAvTMM94lytixcUeSpo4d4Y9/9GTypS/5ML21tXFHJSJZokRSYPbsgUcegU9+0ntvLxht2sCjj8KVV8Ltt8MVV8CBA3FHJSJZoERSYJ57DrZtg899Lu5IGqB5c/j1r+G73/Wn38vKvMJHRAqaEkmBuftuf5K90Ve0H4kZfPObcPPNUF7uO/Pss7rVJVLAlEgKyNtvw9NP+yMZrVrFHU2Gbr8dvvpVH/P93HP9akV1JyIFSYmkgEyf7tMrr4w3jqz5r/+C99+H73zHf779dvjYx2DZsnjjEpHjokRSIHbu9NtaU6d61UJiFBd7Nyq1tZ4pX3sNhg2DUaPgiSdg//64IxSRY1AiKRDf+pZXsn/963FHkiNm8IUvwMqVfqtr3Tq44AI4+WR/YGbmTNi40R+gEZFGxUIT+8ccM2ZMmD9/ftxhHJe1a6F/fxg92gezahLef9/HgH/iCVi69MOf/e//7Z2M9e4N553nnY6JSE6Z2YIQwph6P1Miadz27oUzz/QuUdatgz594o4oBnPmeAbdtAl++MPDP//8533c+PPP97FQRCTrlEhSFFoiuf5675H9wQfh0kvjjqaROHAA3nrLH3D8xjc+POZJr15w6qk+4teQITBunF/OiUhGlEhSFFIi+dWv4Oqr4cYb4c47446mEdu1C7Zvh5/8pP4rlhYtvHnxySfD0KF+9dKjB2zY4ANxnXaaV+p37eo/9+zpn5WXe73MCSdAp07eVHnrVr8sPOkkb4PdurXX74gknBJJikJIJLW1cMcd8J//6fXOM2cWWHcocdu6FR54wOtZVq6EHTugSxe/Nfbaa54ksvW9b9XK62pat/YrIzNvVnfOOb6dCRM8jhC8WXNlpSe9Ll18me7dveO05cs9UXXu7ImsQwfvVmbbNt/GmjW+X8OGQVWVP8hZVuZfjK5dobTU+8/Zts23VVnp8e3b54mvqsqT6bZtXt6pE6xf79vYswc2b/aHRCdPhn//d5g40eMLwT9v2dKXT3XggCfavXs9tooKL+vY0dd/4IDvU02NJ+q9e/1vsGcPtGvnTRFra+HEE2HkSC9fu9bnq672E4CtW/1EoUsXX1+7dvDOO76+zZv9CnT3bh8u9NJLPcbdu/13+vrrvly/fr6u0lLfj8pK/5126eL7UVHhv89du3z5d97x7fbqBQMG+DzV1R7XgQMe8/vv+7q2bfO/84YN/v6EE3w9ffv6Otq1g8GDPd5OnXwft2/3effu9b+fmX9PevTw9+vW+T5UVPjfr6TEr7JXrvTyd97x7/Hq1R5bjx7eQ0RJia/v73/3v+2uXT7kdUmJ/51eegn+/OcGD+egRJKisSeSTZvgppv8OHjZZXDvvf7/LFm0d6//4xcV+T/0++/7P+KyZX7g7NzZD7Bt28Ipp/iBsKLi4JP427fDq6/6gWnHDv+jLV/u5W+95csfbbz6zp39YLt5sx9MwZtB79mTh51PU7Nmhz8c2qePH+h27ID27eG99/y9xKOkxJNIXUKpM2WKf8dfftm/r61b+3d050747/+GG25o0OYSnUjM7DzgLqA58OsQwh1Hm78xJpJdu7wuefp0eOghP+n53vfgllt016RgbdzoB1rwpPT++35m3KbNh+erqfEE1LnzwbP3t9/2f/pmzfzsuazMz36LivwMtVMnb31h5mUbN/rYAief7GfOQ4Z44mvTxg8gvXv7mWzXrgevFrZu9bqjPXv87LxtW786WrfO1/3aa3510KuXf751Kyxa5MsUF3sC7t7df+7a9eAVVps2vmwInojatPF9evddGDTIz8yHDTu4r1u2eEJ9+23/fXTo4Gfxw4f7vB07+rw7d3rSqqry+Pv18+WLiz2+qir/HTdv7vvetq0v17y5/07WrvX1FRX577qiwstat/bfz5Ytvu0OHTyW/fv9QL19ux+M27f3v8HatX5y0bYtrFjh+z9s2MGri7rfa2WlL7Nli+/75s3w4ovet1GPHv5Zhw6+zwsX+naqqnyfBgzwz2tqfN6ZM33/y8q87q9FC6/7a9fu4PeoLuk3y90THYlNJGbWHHgD+DhQDrwK/K8QwoojLZOPRFL3P1R3RV9T4/+vW7f69+mdd/x7vHy5n9guX+7Jo317+Oxn/c5C3RW1iEhjcLREUuh33scBa0II6wDM7AFgKnDERNJQ99wDP/rRwcRwrFc6Onf2MUUuvNCnZ53lJ1siIoWk0BNJb2BDys/lwGmHzmRmVwNXA/Tt27dBG+rSxRv8FBX5q0WLg+/TebVq5Ymje3e/Wu3Rw9epW1ciUugKPZHUdxg+7F5dCGE6MB381lZDNjR1aoMbO4iIJFqh97VVDqQ+610KbIopFhGRJqnQE8mrwAAz62dmLYFpwOMxxyQi0qQU9K2tEEKNmX0J+Cve/PeeEMLymMMSEWlSCjqRAIQQngSejDsOEZGmqtBvbYmISMyUSEREJCNKJCIikhElEhERyUhB97XVEGZWCbzVwMW7AluyGE4h0D43DdrnpiGTfT4xhNCtvg+aXCLJhJnNP1KnZUmlfW4atM9NQ672Wbe2REQkI0okIiKSESWS4zM97gBioH1uGrTPTUNO9ll1JCIikhFdkYiISEaUSEREJCNKJGkys/PMbJWZrTGzm+OOJ9fMrI+ZPW9mK81suZl9Oe6Y8sHMmpvZIjObGXcs+WBmHc3sYTN7Pfpbnx53TLlmZjdG3+llZvZHMyuOO6ZsM7N7zGyzmS1LKetsZs+a2epomrWBvZVI0mBmzYFfAJ8ABgP/y8wGxxtVztUAN4UQPgKMB65tAvsM8GVgZdxB5NFdwNMhhEHACBK+72bWG7geGBNCGIoPPzEt3qhy4rfAeYeU3QzMCiEMAGZFP2eFEkl6xgFrQgjrQgj7gAeARA+8G0KoCCEsjN5X4weY3vFGlVtmVgqcD/w67ljywcxKgI8CvwEIIewLIWyPNaj8KAJam1kR0IYEjqoaQngR2HpI8VRgRvR+BnBRtranRJKe3sCGlJ/LSfhBNZWZlQGjgHkxh5JrPwG+BtTGHEe+nARUAvdGt/N+bWZt4w4ql0IIG4EfAW8DFcB7IYRn4o0qb3qEECrATxSB7tlasRJJeqyesibRbtrM2gGPADeEEHbEHU+umNkFwOYQwoK4Y8mjIuBU4JchhFHA+2TxdkdjFNULTAX6AScAbc3sn+ONqvApkaSnHOiT8nMpCbwcPpSZtcCTyP0hhEfjjifHJgAXmtmb+K3Lj5nZffGGlHPlQHkIoe5K82E8sSTZZGB9CKEyhLAfeBQ4I+aY8uVdM+sFEE03Z2vFSiTpeRUYYGb9zKwlXjn3eMwx5ZSZGX7vfGUI4c6448m1EMItIYTSEEIZ/vf9Wwgh0WeqIYR3gA1mNjAqmgSsiDGkfHgbGG9mbaLv+CQS3sAgxePAFdH7K4DHsrXigh+zPR9CCDVm9iXgr3grj3tCCMtjDivXJgCfA5aa2eKo7D9CCE/GF5LkwHXA/dEJ0jrgX2OOJ6dCCPPM7GFgId4ycREJ7CrFzP4ITAS6mlk58G3gDuAhM7sKT6iXZG176iJFREQyoVtbIiKSESUSERHJiBKJiIhkRIlEREQyokQiIiIZUSIREZGMKJGIiEhGlEhEYmRmXzCzxdGrNuV94nsTkOTQA4kijUA0TsZLIYQT445F5HjpikSkcRgKLI07CJGGUCIRaRyGAcuOOZdII6REItI46IpECpYSiUjjoCsSKViqbBeJmZk1A7YD3UIIe2MOR+S46YpEJH798ZEKlUSkIOmKREREMqIrEhERyYgSiYiIZESJREREMqJEIiIiGVEiERGRjCiRiIhIRpRIREQkI/8fRswQskXf2GgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "tmax = 10\n",
    "tau = 0.5          # transmission rate\n",
    "gamma = 0.05     # recovery rate per node\n",
    "t, S, I = EoN.fast_SIS(BA, tau=tau, gamma=gamma, tmax = tmax)\n",
    "plt.plot(t, S, color = 'r', label='S')\n",
    "plt.plot(t, I, color = 'b', label='I')\n",
    "plt.legend(loc=0)\n",
    "plt.xlabel(\"$T$\")\n",
    "plt.ylabel(\"individuals\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93f9e7ed",
   "metadata": {},
   "source": [
    "#### SIR仿真"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "961fdb5b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Text(0, 0.5, 'individuals')"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYsAAAEICAYAAACuxNj9AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjQuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/MnkTPAAAACXBIWXMAAAsTAAALEwEAmpwYAAA3wUlEQVR4nO3deZzN9f7A8dd7FrMZhrFmMAjJzpC1JlEKIUppX657K9S9/VpvdVtuXbkVbURaVFKSCl0kUhEyluyFbGOJGdswi1k+vz8+Z8xgmHPmnDPfWd7Px+P7+J7zPd/lbbr3vM9nF2MMSiml1PkEOB2AUkqpkk+ThVJKqUJpslBKKVUoTRZKKaUKpclCKaVUoYKcDsBfqlWrZmJjY50OQymlSpWVK1cmGWOqn3m8zCaL2NhYEhISnA5DKaVKFRHZWdBxrYZSSilVKE0WSimlCqXJQimlVKE0WSillCqUJgullFKF0mShlFKqUJoslFJKFUqTxZkmTIDvvnM6CqWUKlHK7KC8Ijl5EsaPh+3bYeNGqFPH6YiUUqpE0JJFfhUqwBdf2KQxciTowlBKKQVosjhbo0bwzDMwYwbMnOl0NEopVSJosijIQw/BRRfBo49CVpbT0SillOM0WRQkKAheegl++w3efdfpaJRSynGaLM6lXz/o3t1WSWVkOB2NUko5SpPFuYjAk0/C/v22/UIppcoxTRbn07MnNGwIkyY5HYlSSjlKk8X5BATATTfBokXw559OR6OUUo7RZFGYIUMgJwemT3c6EqWUcowmi8K0aAHt2tmR3TpITylVTmmyKIwI/PWvsGEDrF/vdDRKKeUITRbuuOYau//2W2fjUEoph2iycEdMjK2Omj3b6UiUUsoRmizcdfXVsGQJHD/udCRKKVXsdIpyd111Ffz3v7Ybbd++TkejVJmXnZNNRnYG6VnppGWmkZ6Vfs4t97xT77MyTvssOycbAIPBGHNqX9Axg+t47rHznOvW9T445snzAebfOp+QoBCf/vfQZOGurl0hLMy2W2iyUOVIjsnh+MnjpGSkcCzjGCkn7f5w2mGS05I5eOIgRzOOnvXlnZaZxonME6RmpnIy+yRZOVlk52Tbvck+9T7/69zPMrMzyczJ9Dr2kMAQQoJCCJRARAQAQRCRU3t3jwmu40U8VhzPDwgIOHXM1zRZuCs0FOLjYd48pyNRyiPZOdkcSjvEgRMHOHDiAH+e+JOk1KSzvvxP7c84fvxk4VWvEcERhAaFEhIUQmhQ6KktIjiCKqFVqBBYgaCAIAIDAu1eztifcTwoIIiw4LBT9wkLCjvtvmHBYYQEnv6s3GfnHq8QWOHUF6/yniYLT1x+OcyZA3v26Cp6qsQwxrA3ZS+r96/m1/2/sv7gev48/icHThzgYOpBklKTyDE5BV4bFBBEpZBKRFaItPuQSKLDo4mNij3r+Jnvq4RWITo8muiwaJ9XeaiSR5OFJ3r3hkcegblz4e67nY5GlROZ2ZkkHktk59Gd7Diyg9+Tf2fHkR3sPLqTfSn72H5k+2nnx0bFUieyDk2im9C1bldqRNQ4baseUZ3q4dWpHFqZkMAQ/fWt3KLJwhMtWkCNGvDdd5oslM+kZaad+vLfeWSn3ed7vTdl72klg0AJpF7lesRGxdKlbheua3YdkRUiiY+Np13tdkSGRDr4r1FllSYLT4jYXlFz59qpP/QXmfJQamYqP+z4gaWJS9l2eBsJexPYkrzlVC8YsFVDdSvVpV7levRo0IP6levbLcruY6NiCQ4MdvBfocojTRae6toVPvoI/vjDrtet1HkcyzjG+gPrWffnOuZsncP8P+aTmpkKQHhwOL0a9uKmFjfRJLrJqYRQu2JtAgMCHY5cqdNpsvBUfLzdz5kDw4c7GooqWbJysthwYAMJexP4cdePLN29lC2Htpz6vF7letzR+g4GXDSAdrXbER0e7WC0SnlGk4WnmjaF5s3tlOWaLMq1pNQk5m6dy+akzfxx+A9+2vUTiccSAYgKjSI+Np7bW99O61qtaRrdlAurXqiNyarU0mRRFH36wJgxkJZmB+qpcmPnkZ3M/n02MzbPYNGOReSYHAIlkPpR9WlTqw3Pxj9L93rdaVS1EQGis+moskOTRVF07w6jR8Py5XnVUqpMys7JZuH2hczbNo9Zv8/i9+TfAWga3ZTHuj7GwGYDaVWzFRUCKzgcqVL+pcmiKLp3tz2hvv9ek0UZZIxh9f7VfLnpSz5a+xE7j+6kQmAF4mPjuS/uPq5oeAUtarRwOkylipVfk4WI/B24BzDAOuBOIBz4DIgFdgA3GGMOu85/HLgbyAZGGmPmuY63Bz4AwoD/AQ8Y4+CydZUrQ/v2sGABPPusY2Eo3zHGsGjHIr7c/CUzf5vJzqM7Abiy0ZWM7jWaPo37EFEhwuEolXKO35KFiNQBRgIXG2PSRGQacCNwMbDAGDNKRB4DHgMeFZGLXZ83By4AvhORJsaYbGA8MAxYhk0WvYE5/ordLb162aqoY8egUiVHQ1FFl5SaxIe/fsiYZWNIPJZIcEAwV114Ff+67F/0adKHGhE1nA5RqRLB39VQQUCYiGRiSxR7gceBeNfnk4FFwKNAf+BTY0wGsF1EtgIdRWQHUMkYsxRARD4EBuB0soiPh//8x7Zb9OrlaCjKM8dPHuez9Z/x1W9fMXfrXLJysris/mU8F/8cN7a4kbBg7bSg1Jn8liyMMXtE5GVgF5AGfGuM+VZEahpj9rnO2SciuT/d6mBLDrkSXccyXa/PPH4WERmGLYFQr149X/5zztaxo223WLJEk0UpcST9CK/8/AqjlowiKyeL2KhYhncYzu1tbqdNrTZOh6dUiebPaqgq2NJCA+AI8LmI3HK+Swo4Zs5z/OyDxkwEJgLExcX5t00jKgri4uyU5c8849dHKe/sPLKT/yz+Dx+v/ZjUzFT6X9SfER1HcHns5TruQSk3+bMaqiew3RhzEEBEZgBdgD9FpLarVFEbOOA6PxGom+/6GGy1VaLr9ZnHndejB7z6KmRkQIhO0VyS5Jgcpm+czthlY1mWuIyggCAGXzyYhzo/RPsL2jsdnlKljj9HDe0COolIuNifb1cAm4CZwO2uc24Hvna9ngncKCIhItIAaAz84qqyShGRTq773JbvGme1bw+ZmbB2rdORqHy2JG+h06RODJk+hOS0ZJ6Nf5ZtI7fxyaBPNFEoVUT+bLNYLiLTgVVAFrAaW0VUEZgmIndjE8r1rvM3uHpMbXSdf7+rJxTAveR1nZ2D043buTp1svtly6BDB2djUSSnJvP0908zcdVEIoIj+HDAh9zc6mYdSa2UD4iTwxX8KS4uziQkJPj3IcZA/fq2sXv6dP8+S51Tdk42b614i2cWPcPh9MPc2eZOXujxArUjazsdmlKljoisNMbEnXlcR3B7QwSuuAJmzdL1LRyy9s+13DPzHlbsXUHPhj0Z3XM0bWu3dTospcocLZ97q3t3SE6GzZudjqRc+T35d26ZcQvtJrRjx5EdTB00lW9v+VYThVJ+oiULb116qd0vWADNmjkbSzlwLOMYTyx4ggkrJxAgAYy8ZCT/7P5PXRtCKT/TZOGtRo2gSRP45htd38KPjDFM2zCN4XOGk5yazN/i/sZTlz6l7RJKFRNNFt4SseMtpkyB7GwI1OUwfW3boW0Mmz2MhdsX0qZWG+bcPIe4C85qf1NK+ZG2WfhC166QkgIbNjgdSZmSkpHCHV/dwUVvXcSSXUt45cpXWHb3Mk0USjlASxa+0KWL3S9dCq1aORtLGbHp4CYGfDaArYe2MrLjSB7p+ohWOSnlIC1Z+EKDBlCrFvz4o9ORlAmfrf+MDu904HDaYRbctoAxvcdoolDKYVqy8AUR2ytKk4VXjqYf5e/z/s77a96nS90uTBs8jTqVCpxgWClVzLRk4Stdu0JiIuze7XQkpdLhtMMMmT6E99e8z0OdH2LR7Ys0UShVgmjJwlc6d7b7pUuhbt3zn6tOczjtMD0/6snqfauZ0HcCw9oPczokpdQZtGThK61bQ2ioTRbKbdsPb+fyyZez7s91zLppliYKpUooTRa+UqGCXQxJk4XbZv8+m27vd2Pn0Z3MumkWfZr0cTokpdQ5aLLwpc6dYdUqSE93OpISb0LCBPpN7Ud0WDQLb1vIVRde5XRISqnz0GThS1272sWQVq50OpISbcamGdz7zb30adyHX/7yi07+p1QpoMnCl/I3cqsC/bz7Z26ecTOXxFzCtOunERoU6nRISik3aLLwpRo1oGFDu3KeOsu3277l8smXUyeyDl8N+Yrw4HCnQ1JKuUmTha916qQliwJ8vuFzrp5yNU2im7DsnmXUrFjT6ZCUUh7QZOFrHTvC3r2wf7/TkZQYMzbNYMj0IVxS5xJ+vutnqoVXczokpZSHNFn4Wrt2dv/LL87GUUJMXTeVGz6/gbgL4ph/63wiQyKdDkkpVQSaLHwtLg6Cg2HJEqcjcdy4FeMYOmMonWI6seC2BURUiHA6JKVUEWmy8LWwMFu6+PlnpyNx1KzfZjFizgj6NenHwtsXaolCqVJOk4U/dO0KK1ZARobTkThi1b5V3PjFjbSr3Y6pg6ZSIbCC0yEppbykycIfuna1iWL1aqcjKXa7j+6m39R+VAuvxqybZmnVk1JlhCYLf+ja1e4XL3Y2jmKWmZ3J4M8Hk5KRwuybZlOrYi2nQ1JK+YhOUe4PNWva1fOWL3c6kmL15MIn+WXPL3x+/ee0rNnS6XCUUj6kJQt/KWeD877d9i2jfx7NX9v/lcEXD3Y6HKWUj2my8JcuXWDPHti1y+lI/G7/8f3c+uWtNK/enDFXjXE6HKWUH2iy8JfcwXlr1zobh5+lZ6Vz/efXcyzjGJ8N/oyw4DCnQ1JK+YEmC39p2RJEICHB6Uj8JsfkcPtXt7N412LevfZdmtdo7nRISik/0WThL5GR0KxZmV7b4u2Et5m2YRov9XyJoS2HOh2OUsqPNFn4U1yc7RFljNOR+NyKPSt46NuH6NmwJw93edjpcJRSfqbJwp+6dIGDB+GPP5yOxKcOpR1iyPQh1IyoyZTrpiAiToeklPIzvyYLEYkSkekisllENolIZxGpKiLzRWSLa18l3/mPi8hWEflNRK7Kd7y9iKxzffa6lJZvp8sus/vvvnM2Dh/Kzsnm+s+vZ0/KHj4Z9Ak1Imo4HZJSqhj4u2TxGjDXGHMR0BrYBDwGLDDGNAYWuN4jIhcDNwLNgd7AOBEJdN1nPDAMaOzaevs5bt9o2hTq1oVvv3U6Ep8ZtXgUC7cvZNw14+hSt4vT4SilionfkoWIVAIuBd4FMMacNMYcAfoDk12nTQYGuF73Bz41xmQYY7YDW4GOIlIbqGSMWWqMMcCH+a4p2URs6aKMjOTenLSZZ354hhtb3Mhdbe9yOhylVDHyZ8miIXAQeF9EVovIJBGJAGoaY/YBuPa59Rh1gN35rk90Havjen3m8bOIyDARSRCRhIMHD/r2X1NU7dvbwXl79zodiVfSs9IZMn0IUaFRjLlqjLZTKFXO+DNZBAHtgPHGmLbACVxVTudQ0LePOc/xsw8aM9EYE2eMiatevbqn8fpH5852X8qn/hi3Yhxr/1zL5AGTdYJApcohtyYSFJHrsW0PKSLyJDYJ/NsYs+o8lyUCicaY3DqY6dhk8aeI1DbG7HNVMR3Id37dfNfHAHtdx2MKOF46tG0LISF2MaRBg5yOpkjW7F/D098/zeWxl3NN42ucDkepEi0zM5PExETS09OdDuW8QkNDiYmJITg42K3z3Z119iljzOci0g24CngZ2+h8ybkuMMbsF5HdItLUGPMbcAWw0bXdDoxy7b92XTIT+EREXgUuwDZk/2KMyRaRFBHpBCwHbgPecDNu51WoAB06lNqV805mn+SmL26icmhlplw3xelwlCrxEhMTiYyMJDY2tsRW1xpjSE5OJjExkQYNGrh1jbvVUNmufR9stdLXgDvLn40ApojIWqAN8CI2SfQSkS1AL9d7jDEbgGnYZDIXuN8Yk/vce4FJ2EbvbcAcN+MuGbp0sSO5S/gvjYK8/PPLbE7azIS+E6gdWdvpcJQq8dLT04mOji6xiQJARIiOjvao9ONuyWKPiEwAegIviUgIbiQaY8waIK6Aj644x/kvAC8UcDwBaOFmrCVPly4werSdJ6pbN6ejcdsve37hqe+fYvDFg+nbpK/T4ShVapTkRJHL0xjdLVncAMwDeru6v1YFdI4Hd3XvDgEBsHCh05G4LSMrg7u+vovaFWszqd8kp8NRSjnsvMnCNdq6KhAKLAKSXe8zgLI7naqvVa0KrVqVqpHcD337EBsObmBC3wlUDq3sdDhKKQ+88MILNG/enFatWtGmTRuW+2CsV2HVUCvJ6756ZjdWgx1LodzRuze8/DKkpkJ4uNPRnNeSXUt4a8Vb3Bd3H32a9HE6HKWUB5YuXcrs2bNZtWoVISEhJCUlcfLkSa/ve96ShTGmgTGm4Rn73E0ThScuvRSyskp8r6jjJ49z59d3Ur9yfV7q9ZLT4SilPLRv3z6qVatGSEgIANWqVeOCCy7w+r7uNnDjmvCvMbZKCgBjzI9eR1BedOsGQUGwYAH07Ol0NOf07KJn2XJoC9/f/j0VK1R0OhylSrcHH4Q1a3x7zzZtYOzYc3585ZVX8txzz9GkSRN69uzJkCFDuCx3UlMvuNXALSL3AD9iG7mfde2f8frp5UlkpJ36Y8kSpyM5p/UH1vPqsle5p+09xMfGOx2OUqoIKlasyMqVK5k4cSLVq1dnyJAhfPDBB17f192SxQNAB2CZMeZyEbkImzSUJzp2hPfeg5Mn7WC9Eubh+Q9TKaQSo3qOcjoUpcqG85QA/CkwMJD4+Hji4+Np2bIlkydP5o477vDqnu52nU03xqQDiEiIMWYz0NSrJ5dHPXrAiRMlcp6oH3f+yNytc3my+5NEh0c7HY5Sqoh+++03tmzZcur9mjVrqF+/vtf3dbdkkSgiUcBXwHwROUxpmp+ppIiPt+Mt5s3LWxipBDDG8NT3T1EltAr3drjX6XCUUl44fvw4I0aM4MiRIwQFBXHhhRcyceJEr+/rVrIwxgx0vXxGRL4HKmOn5FCeiIqySeKbb+DFF52O5pTxCeP5ceePjLtmHOHBJbtbr1Lq/Nq3b8/Pfuh16W4Dd73cDdgOrAF0nuqi6NED1q6F5GSnIwHsetr/XPhP4mPj+WvcX50ORylVQrlbDfUNeYPyQoEGwG/YJVCVJ+Lj7f6nn2DAACcjAeCR+Y9wNP0or/V+jQDx9yq7SqnSyq1vB2NMS2NMK9e+MdARWOzf0MqoDh0gLKxEzBO1/sB63l39Ln/v9Hda1WzldDhKqRKsSD8lXYsedfBxLOVDSIgdzT1vnqNhGGN4YO4DVAmtwhPdn3A0FqVUyefuSnn/yPc2ALtSXglZ5LoU6tcPhg+HzZvhooscCeGjtR+xcPtCxl0zTrvKKqUK5W7JIjLfFoJtw+jvr6DKvGuvtfuZMx15fGZ2Js//+Dxta7XVRm2llFvc7Tqro7V9qW5daNkS5syBRx4p9sd/sOYDth7ayldDvtJGbaXKoIoVK3L8+HGf3vO8yUJEZmF7QRXIGHOtT6MpT3r3tlMBnDgBERHF9tiMrAye//F5OsV04tqm+p9PKeWewn5Wvgy8gh1bkQa849qOA+v9G1oZd/XVkJlpSxfF6J1V77D72G7+ffm/S8XSj0qpkuG8JQtjzA8AIvK8MebSfB/NEhGdntwb3btDtWrwxRcweHCxPDI9K50XfnqB+Nh4ejToUSzPVKo8c2CGcr9xt8K6uoicWuxIRBoA1f0TUjkRFAQDB8Ls2ZCeXiyP/Hjtx+w/vp8nuz+ppQqllEfcHcH9d2CRiPzheh8LaDcabw0aBO+8YxdE6uPf5UtzTA6vLn2VtrXaaqlCqWLi0AzlfuFub6i5ItIYyB0UsNkYk+G/sMqJ+Hi7KNLMmX5PFnO2zGFT0iY+HvixliqUUh47bzWUiPRw7a8D+gCNXFsf1zHljZAQ29D95ZeQne23x5w4eYLHFjxGTKUYbmh+g9+eo5QquworWVwGLAT6FfCZAWb4PKLyZvBgmDYNFi/22xoXry1/jfUH1jP9+ukEBwb75RlKqZLD12MsoPDeUP9yvbzHGOO/n77l2dVX2xLGl1/6JVkcyzjGf3/+L32b9GXQxYN8fn+lVPngbm+o7SIyUUSuEK3w9q2KFeHKK2HGDMjJ8fntX1v2GkfSj/DMZc/4/N5KqfLD3WTRFPgOuB+bON4UkW7+C6ucueEG2L3bVkX50OG0w7y67FX6NelH+wva+/TeSqnyxd31LNKMMdOMMdcBbYFKwA9+jaw8GTgQwsPhk098ettXlr7C0fSjPHf5cz69r1Kq/HF7FjkRuUxExgGrsKvlabcaX4mIsG0Xs2eDOedUXB45nHaY15e/zuCLB9OmVhuf3FMpVX65uwb3duBB4CeghTHmBmPMF/4MrNy59lrYswd8tND62GVjSTmZwlOXPuWT+ymlyjd3SxatjTEDjTFTjTEn/BpReXXddbaE8cEHXt/qwIkDvLz0ZQY1G0TLmi29j00pVaoEBgbSpk0bWrRoQb9+/Thy5IjX9yxsUF7uYgv/FpHXz9y8frrKk9srau5cr3tFvZ3wNqmZqTx/+fM+Ck4pVZqEhYWxZs0a1q9fT9WqVXnrrbe8vmdhJYtNrv3Kc2yFEpFAEVktIrNd76uKyHwR2eLaV8l37uMislVEfhORq/Idby8i61yfvV5mu+8OGgSJiV5VRaVlpvHmL2/Sp3EfmlVv5sPglFKlUefOndmzZ4/X9ylsUN4s136yF894AJt0KrnePwYsMMaMEpHHXO8fFZGLgRuB5sAFwHci0sQ1GHA8MAxYBvwP6A0U70IQxeHaayEsDKZOhW5F65k8+dfJHEw9yMNdHvZxcEopTz0490HW7F/j03u2qdWGsb3HunVudnY2CxYs4O677/b6uYVVQ80SkZnn2gq7uYjEYOeUmpTvcH8gN/lMBgbkO/6pMSbDGLMd2Ap0FJHaQCVjzFJjjAE+zHdN2RIZaRPGtGl2YSQP5ZgcXl/+Ou1qt+PS+pcWfoFSqkxKS0ujTZs2REdHc+jQIXr16uX1PQubG+pl1/46oBbwsev9TcAON+4/FngEiMx3rKYxZh+AMWafiNRwHa+DLTnkSnQdy3S9PvP4WURkGLYEQr169dwIrwQaOhQ++wzmzYO+fT269KvNX7EpaRNTB03VmWWVKgHcLQH4Wm6bxdGjR+nbty9vvfUWI0eO9Oqe5y1ZGGN+cK2W19YYM8QYM8u1DQXOW08iIn2BA8YYt9o2gIK+3cx5jhcU70RjTJwxJq569VK6NtPVV0ONGvDeex5fOj5hPDGVYrj+4uv9EJhSqrSpXLkyr7/+Oi+//DKZRaityM+fK+V1Ba4VkR3Ap0APEfkY+NNVtYRrf8B1fiJQN9/1McBe1/GYAo6XTcHBcNNN8M03kJLi9mUbDmzguz++4764+wgMCPRjgEqp0qRt27a0bt2aTz/91Kv7uJssclfKWyQii4DvsYP0zskY87gxJsYYE4ttuF5ojLkFmAnc7jrtduBr1+uZwI0iEuJKRo2BX1xVViki0snVC+q2fNeUSdlDhvL+yaFkv/MeS5a415N20qpJBAUEcXc77xuylFKl25lTlM+aNYtbb73Vq3s6sVLeKGCaiNwN7AKudz1jg4hMAzYCWcD9+aZFvxf4AAjD9oIqez2hXGxTQ0egI3c9dPbnS5ZAly6nHzuWcYz317zP4IsHUyOixtkXKaWUl9xdgxugPXbt7SCgtYhgjPnQnQuNMYuARa7XycAV5zjvBeCFAo4nAC08iLXUWbYMOncu/LyuXe0+MhIOHYKgIDsI72jGUf6v8//5N0ilVLnl7txQH2F7RnUDOri2OD/GVW5s3mxLE/kTxeHDkJOVw5Sa/yA77hKMgdGjT78uJcU2b2z8PZ0xy8bQq2EvnYZcqRLC+GhCUH/yNEZ32yzigK7GmPuMMSNcm3f9sBSLF0OzMwZZGwNRUSCBAQz9ZwMCEn6BhAQefth+ZgycPJl3fvObP2D/8f083u3xYo1dKVWw0NBQkpOTS3TCMMaQnJxMaGio29eIO/8gEfkcGJk7PqI0iIuLMwkJCU6HcU5r10Lr1nnvjx2zVUunOXoU6tSx04BMPnsQvYiBe1tDdgU6rVvB0p91bIVSTsvMzCQxMZH09HSnQzmv0NBQYmJiCA4OPu24iKw0xpxVc+Rum0U1YKOI/AKcatg2xlzrTbDlVWbm6YninPm6cmW48054+2146SWoVeu0j+dtnc9VH6+DWRNYtlJ49FF7mlLKOcHBwTRo0MDpMHzO3WqoZ7BTbLwIvJJvU0UQly9nF1qwGzECsrLg9bMn+f3P4v9QJ7IO++feAdh2jb/8xXdxKqVULneXVf2hoM3fwZVFa9fC+vX2dXb2+c8FoEkTuP56ePNN2/LtsuHABhbtWMT9He6nZrUKvP++PT5pEvhgzjCllDpNYRMJLnbtU0TkWL4tRUSOFU+IZUdWFtxwA0RHQ1ISBLhbrnviCdv96bXXTh167sfniAiO4C/tbVHijjsgd8r6994r0mwhSil1ToXNDdXNtY80xlTKt0UaYyqd71p1ti5d4LffYMwYmzDc1qYNDBgAY8fC0aPsOrqL6Runc2/cvVQLr3bqtPvugxUr7Ou777ZdcktwhwylVCni7m9b5aWMjLwv8qFDi3CDp56yvaPGjmX8ivEADO84/KzT4uJsM0eugABNGEop72myKCa53ZmffDJ3Sg8PtWsHffpwfOKbTFg5gf5N+1M/qn6Bp77++untIffcU4TnKaVUPposikH+tdKfe86LGz30EO/WTeJw+uFCV8ILCMhbP+m992yC2rXLi2crpco1TRbF4I478l57syZR9qXdef2yEDrvD6Jz1VaFnh8UdPqCe/Xrnz76Wyml3KXJohh87ZpQ3Z2pxs9n5u+z+CM8g4d+yOJUX9lCBAXZXli54nRGL6VUEWiy8LPt2+1+0CDvShUAry57ldioWAZEd4XHH4c9e9y6LjAwL1GtWwc//uhdHEqp8keThZ81dK0vOGqUd/dZsWcFi3ct5oFLHiDwg8m2uDBihNtdnUTgwAG7v+wyuymllLs0WfjRkiV5ry+80Lt7jVk2hsgKkdzV9i5o1Mi2lH/5Jbzxhtv3qF4dFi60r3/8EV44a+UQpZQqmCYLP+rWze69XPqWgycO8vnGz7mn3T1UCnGNhXzoIejXz+6XLnX7XvHxeY3cTz4Jw88eqqGUUmfRZOEn+WcnHjLEu3t9tuEzsnKyuKPNHXkHAwLstOV169o5RA4edPt+wcF2JDnYKUK0DUMpVRhNFn4yx7VK+Nix3t3HGMO4FeNoV7sdLWu0PP3DKlXgiy9sohg61M2ZCa0mTfJGlF92Gezf712cSqmyTZOFnyxaZHsheTt6etGORWxK2sSIjiOQgrpTtW1rZ6T97jt48UWP7h0XB7Nm2Ubv3r1PHzyolFL5abLwg5QUO+XGwIEQEeHdvd5e+TZRoVEMaX6euqy777ZVUc8/DwsWeHT/vn1tKejXX21BJTXVu3iVUmWTJgs/6N/f7nv18u4++1L28dXmr7it1W2EBYed+0QRGD8eGjSwRZnkZI+ec9VVcMst9rW3yU0pVTZpsvCD77+3e29XrRufMJ6T2Se5v+P9hZ9ctaqdBGrvXlvS8HC4+EcfQbNm9vU77xQhWKVUmabJwsfS0uy+WTPvRmwnpSbx8s8v07dJX5pEN3Hvoq5dbbvF11/DP//p8TMTEux+2LC8hZSUUgo0WfjczJl27+2X7eQ1k0nLSmPUFR4O/f7HP+zMhaNGedx+ER5uO1YFBNjxF48+6tmjlVJllyYLH3vlFbu/9FLv7vPVb1/RokYLmtdo7tmFIjZTNW1qB3isXOnR5dWq2UbuWrVg9Gi44ALPHq+UKps0WfhQdjasXg2RkbbbbFFtOLCBxbsWc0vLW4p2g/BwmD0bQkJsK3vugAo3hYTAjh329b59Nv94O2OuUqp002ThQ23b2vn93nvPu/uMWjKKsKCw00dse+rCC2HxYoiKgh497DgMD4SEnL4WRmAgJCUVPRylVOmmycKH1q2z+9yus0Wx48gOPln3CffG3UvNijW9C6hBA5sw6teHa689fWZDNwQF2RJFbimpenU4ccK7kJRSpZMmCx/J/RKNj7dzLxXVxJUTEYQHOz3oi7Bso8P8+VCjhp3XY/x4jy4XsaWljh3t+4oVdbU9pcojTRY+kvsd/PTTRb9HZnYmk3+dzBUNr6Bu5bq+CQygdm1YswYuvxzuuw+eesrtdTByLV8OjzxiX3fr5tE0VEqpMkCThY/krnLqTS+oqeunsjdlLw9e8qBPYjpNVJSd1+POO+Hf/7brYXjopZfy2svj4uD4cd+HqZQqmfyWLESkroh8LyKbRGSDiDzgOl5VROaLyBbXvkq+ax4Xka0i8puIXJXveHsRWef67HUpcEY952RmwsaNcNtt3vWCGrdiHE2jm9L7wt6+Cy6/oCCYNAkGD4ZnnoH77/doanOAb7+1JYs1a+zs6G6u7KqUKuX8WbLIAh4yxjQDOgH3i8jFwGPAAmNMY2CB6z2uz24EmgO9gXEikvvVOx4YBjR2bX76Ni2a5q6hELndTYti6e6lLN+znOEdhxc8u6yvBATAhx/CyJHw9ts2+Lfftg0TbvrpJ5g61c5S26GDx/lGKVUK+S1ZGGP2GWNWuV6nAJuAOkB/YLLrtMnAANfr/sCnxpgMY8x2YCvQUURqA5WMMUuNMQb4MN81JcKWLXY/aVLR7/Hvn/5N9fDq3N76dt8EdT5hYfDaa3ZQSP36cO+9tgX7p5/cvsWNN8KMGXbOwkaNbBu6UqrsKpY2CxGJBdoCy4Gaxph9YBMKUMN1Wh1gd77LEl3H6rhen3m8oOcME5EEEUk4WEw/d3MneL3ySmjcuGj3WH9gPf/b8j/u73A/kSGRvguuMK1a2QaITz6BP/+0DS4DB+Yto1eIgQPh558hI8N2F85NmkqpssfvyUJEKgJfAA8aY46d79QCjpnzHD/7oDETjTFxxpi46tWrex5sEUyfbvcvvVT0e0xaNYnggGD3Zpf1h5tust/0L7xg55Nq3tyWNtyoV2vf3i70lJZml9RISfF7tEopB/g1WYhIMDZRTDHGzHAd/tNVtYRrf8B1PBHI3180BtjrOh5TwPESYcoUuPhiaN26aNfnmBw+WfcJ/S/qT7Xwar4NzhPh4fDEE7B1q00U77xji0oPPADHzpfjoXNn+OorOyixRw8d6a1UWeTP3lACvAtsMsa8mu+jmUBuxfztwNf5jt8oIiEi0gDbkP2Lq6oqRUQ6ue55W75rHLVxo63mHzSo6NORL961mIOpBxl40UDfBldUNWrAG2/AH3/Y2WvfeAMuusgu3XqeRvD+/eHLL2H9eujeHXbvPuepSqlSyJ8li67ArUAPEVnj2q4BRgG9RGQL0Mv1HmPMBmAasBGYC9xvjMkd+nUvMAnb6L0NmOPHuN22aJHdX3550e8xZtkYosOi6d/UizlC/KFePVu6WLIEGjaEESPsTLZvvnnOkka/fjBvnl1/qVs3+P33Yo5ZKeU3YjwcyVtaxMXFmYTc1Xz85MILYedOO/1FUUoW2w5to8mbTXi066O8eMWLvg/QV4yxs9j+61+2B1VICFx/vV07o02bs/7xq1fbpVrBJo+2bYs/ZKVU0YjISmNM3JnHdQR3ESUlwbZt0Lt30aug/rXoX4QFhTG843DfBudrIrbYsGoVLFtmF+z+4gto1862a7z9th2d55rHvG1b+OEHqFDBntKzp+0xpZQqvTRZFFHuOttPPFG065NSk5i+cTp3trmTCyJL0QpDl1xiB5Ts2mXbMyIjbYN4TIydlvb66+GDD2gWup0lS2xvqQUL7PQgWi2lVOmlyaKIvvzSrirXoUPRrn/zlzfJyM7g3g73+jaw4lKtml17ddUq+OUX25bRr5+dEv3OO6FhQ+oPaEtCz8f44NFN7N1raNcOPvrI6cCVUkWhbRZFkJZmf0QPHQoTJ3p+/YmTJ6g7pi7d6nVj5k0zfR+gk7KzYdMmmDXLbitWQFYWicENGBr6BT+ltOXaDvuYOqsi4TWLcQCiUsot52qzCHIimNJuwgS7fsWgQUW7/vONn3M4/TAPd3nYt4GVBIGB0KKF3R5/HA4fhp9+Iubnn1n404OMXHYT41f8jaq10tnfKI6oRtG2a26rVnbftq0d86GUKlG0ZFEEuQ3aJ096vtCRMYb2E9uTkZ3B+nvX+3fSwJIoNZVbBqQwZb5dBXBVs5tpu+NLW1zLFRlp1+CIibFzVrVoAbVq2aqvBg2gUiWHgleq7NOShY9s3Zr3uigr4i3YvoDV+1czoe+E8pcoAMLD+fjbcIb+D/r0gfabp/DM0zk8fMNOwn7/FTZsgAMH7GCN7dvhv/89e6WlihXt/OjR0dCkCdSpA1Wr2gGFNWvapBIWZrtjVagAoaF2PY8AbaJTqqg0WXgody6oohZaxiwbQ9WwqsUzu2wJds01du7CYcPgX88G8NnnDfjmmwbEDhhw+okZGTZDJyfbC7Zvtz2x9u2z/Ze/+cYml8JKyCI2oVSqZBNJZKRNLpUr2+NRUTbJ1KiRd7xiRVsllpt0crcg/b+NKn/0f/UemjnTVqu3b+/5tfO3zed/W/7H6J6jCQkK8X1wpUyNGrZX2ezZcPPN9u/66qu2M9UpISF5C4acS06OXVzjzz9t4khKsknm5Em7pabCoUP289RUW+V19Kidk2TtWtuuUsj8V6cJDLSllbAwuw8KyivBhIbaJBMZaT8Xsa8rVrT/lgoV7D442F6Xuw8KssfDw20JKCDAXpv72p33FSrY2PIfO9eWe15wsL2uoFLuuUq+3p7ri3uUhHPLGU0WHti9G5Yuheef9/xaYwyjloyiZkRNRl4y0vfBlVK54/0WL7bLg991F3z6qe1lVr++mzcJCLClg6pVoVmzogWSnW1LLwcO2NWcjh2z68aeOGGXQsxNPBkZdktPt0knPd3OmZX7WXq6vTYpySamAwfsF3Pu9ZmZRYtPlR4lIZEdPWp/uPiQJgsPvPKK3Q8swpx/CXsTWLh9oZYqzqFVK/juOztL+ssvQ2wsPPggjBlTTAEEBuZVQflTTo5NGllZNnHk3+cmIGPsebnb+d7nvs7OtvfNf15BW3b26a9zE+GZzlWtV9Bxf51b3M8rS7F5s77zOWhvKA/kJvCi/Mmu/OhKVu5bybaR24gKjfJpXGXNrl22mi8pyTaCP/AAXHGFtk8rVRx0bigvbdtm98OGeX7twu0Lmf/HfB7v9rgmCjfUqwf798OLL9qV+K680iaPyZNPTT+llCpmmizcNMO1dNPjj3t2XXZONn+b/TcaRDXgvg73+T6wMiow0P6td+6EceNsKeOOO+xCU19/XbTSnVKq6DRZuOmLL+yv29hYz66bsWkGWw5t4b+9/kt4sI5M9lTuPIW//QZjx9pOTwMG2NlsP/3UtiErpfxPk4UbEhNh+XK47jrPrjPG8NKSl7iw6oUMuGiAX2IrL8LDbdvF7t12pt/UVLt0eESEnQJdk4ZS/qXJwg1vvmn3niaLBdsXsHLfSh7p8giBAb7vnVAeBQfbHlMbN9reU3Fxdgr0iAg7Y3pKitMRKlU2abIoRFoavPSSfX3RRZ5dO2bZGGpG1OS21rf5PrByLjDQ9pBasQJ+/NF2vR050s788fzzdqiEUsp3NFkUYto0u//uO8+uyx2tfX+H+3VchZ917w5r1tieUx06wNNP27F506drQ7hSvqLjLM4jJwdatrS/Ytescb+ff0ZWBi3GtyDH5LD2b2uJqBDhVRzKM2vXwq232n3TptC3r10Jtk0bpyNTquTTcRZFMHu2rRt/7DHPBoS9t/o9th7aymu9X9NE4YBWrWDlSrsqX3S0HQWeO5/XhAl2sLRSyjNasjiHjAzbaJo7i4K705H/nvw77Sa0I+6COBbevpAA0XzstN27YfRou3Dfzp12aYw777Sjwzt2LNpU80qVVVqy8NDf/24TxRtvuP9lkp2TzV1f30VwYDAfX/exJooSom5d+99x+3Y7o3mHDrbTQrdudnnca6+F+fPPXjZDKZVHJxIswNixMH687W0zfLj7183YNIMlu5cwqd8kYirF+C0+VTQidh2Na66xE8z+8INNHrNn21KHiB2z0bkzXHYZdOni84k7lSq1tBrqDEuW2F+ctWrZNXci3GxyOJx2mObjmlM1rCpr/raGoADNw6VFejq8/TZ8/z388Ydtp8rJsQMBe/Sw1VUDB9pF+JQq685VDaXJIp/MTDudx9699kujQQP3rkvPSufqKVezeNdilt29jPYXFGFlJFViHDtmx27Mmwdz5thJJEVsSeOqq2yJs1Uru56RUmWNJgs3LV4MtWtDo0buX/Pg3Ad5bflrvHvtu9zV9i6Pn6lKLmPssuDTp9sqq5Ur88ZuXHONnRG3dWv7I6NePZ1GXZV+miz8ZOq6qQydMZThHYbzxjVv+P15yllJSfDTT/DWW3bdjS1b8j6rXNlWWzVpYgcFtm9vR/3rkt2qNNFk4Qfzt83nmk+uIe6COL6//XtCg7Q1tLzZvRt+/91WWy5bZhvNd+3KWz01LMxWWTVpApdcAp062ffaXVeVVJosfGz6xuncPONmGldtzE93/kSVsCp+e5YqXbKybOeIlSvt3FXr1tlG8/377edhYbbU0batbQdp1sxWe2obiCoJNFn40OvLX+eBuQ9wSZ1LmHXTLKpHVPfLc1TZYYwtcSxblrf9+qudqDJX1aq27aNuXTshYv36dqtZ0x6LiYEQnWZM+dm5koXWpnrgZPZJRi8ZzVPfP0W/Jv2Ydv00rXpSbhHJ+/IfMsQey8qC9ettKWTLFlultX277X31ww92oaczRUfbRNKoke3eXaOG3WJi7LGGDW3JRSlf02ThpmWJyxgxZwQJexMYfPFgPrnuE4IDteJZFV1QkJ3c8FwTHB47ZqcnOXDAJpLdu2HfPltC2bzZdu9NTj77ugsusIkjNtYml+hom0AaNbIj1qtVs1uVKtp7S7mv1CQLEekNvAYEApOMMaP8/UxjDKv3r2b0ktF8tuEzaleszccDP2Zoy6GIiL8fr8q5SpXsrMfnk5Vl1+7YtcuWSPJvP/xgk8mJEwVfGxBgE0Z4uN0iIjzboqJsIqpe3b4ODdWeX2VZqfhPKyKBwFtALyARWCEiM40xG/3xvJ1HdjJl3RQ+Xf8p6w6sI7JCJP/X+f948tInqRxa2R+PVKpIgoLsuKDatW1vq4JkZNjG9eRk2/X34MG8/aFDdgR7aqpNKidO2GO7d+e9P3789LaV8wkMtEkjJMSWZnKTUIUK9rPcLSjo9PcFfRYcnJfIcu91vi3/MwICTt8XdkzEbpD3Ov8xVUqSBdAR2GqM+QNARD4F+gM+TRbGGPpN7cc3W76xD63TkbeueYuhLYcSFRrly0cpVWxCQvLaS4oqJ+f0hHLihG1TSU621WTHjtmkk39LS8s7NzPTTtSYu6Wlnf4+/5aVZfeZmfa81FT3k5U/FZREzjzm6XF/nbtune/nNSstyaIOsDvf+0TgrN9RIjIMGAZQr149jx8iIjSu2pinL32aO9veSWxUbNGiVaqMCQiwXXud6t6bk5NXAsq/5SaT1FRbgsrOtufm35/rdf5juZ1CjTl78+R4STnXH21RpSVZFFQYPKvPrzFmIjARbNfZojxoTO8xRblMKeVHAQF51U3KGaWlL0QiUDff+xhgr0OxKKVUuVNaksUKoLGINBCRCsCNwEyHY1JKqXKjVFRDGWOyRGQ4MA/bdfY9Y8wGh8NSSqlyo1QkCwBjzP+A/zkdh1JKlUelpRpKKaWUgzRZKKWUKpQmC6WUUoXSZKGUUqpQZXY9CxE5COws4uXVgCQfhlOa6d/C0r9DHv1b5CmLf4v6xpizFukps8nCGyKSUNDiH+WR/i0s/Tvk0b9FnvL0t9BqKKWUUoXSZKGUUqpQmiwKNtHpAEoQ/VtY+nfIo3+LPOXmb6FtFkoppQqlJQullFKF0mShlFKqUJos8hGR3iLym4hsFZHHnI7HKSJSV0S+F5FNIrJBRB5wOianiUigiKwWkdlOx+IkEYkSkekistn1v4/OTsfkBBH5u+v/G+tFZKqI+HgR05JHk4WLiAQCbwFXAxcDN4nIxc5G5Zgs4CFjTDOgE3B/Of5b5HoA2OR0ECXAa8BcY8xFQGvK4d9EROoAI4E4Y0wL7LIJNzoblf9pssjTEdhqjPnDGHMS+BTo73BMjjDG7DPGrHK9TsF+IdRxNirniEgM0AeY5HQsThKRSsClwLsAxpiTxpgjjgblnCAgTESCgHDKwcqdmizy1AF253ufSDn+gswlIrFAW2C5w6E4aSzwCJDjcBxOawgcBN53VclNEpEIp4MqbsaYPcDLwC5gH3DUGPOts1H5nyaLPFLAsXLdr1hEKgJfAA8aY445HY8TRKQvcMAYs9LpWEqAIKAdMN4Y0xY4AZS7tj0RqYKtdWgAXABEiMgtzkblf5os8iQCdfO9j6EcFC3PRUSCsYliijFmhtPxOKgrcK2I7MBWTfYQkY+dDckxiUCiMSa3lDkdmzzKm57AdmPMQWNMJjAD6OJwTH6nySLPCqCxiDQQkQrYBquZDsfkCBERbL30JmPMq07H4yRjzOPGmBhjTCz2fxMLjTFl/ldkQYwx+4HdItLUdegKYKODITllF9BJRMJd/1+5gnLQ0F9q1uD2N2NMlogMB+Zheze8Z4zZ4HBYTukK3AqsE5E1rmNPuNZBV+XbCGCK6wfVH8CdDsdT7Iwxy0VkOrAK23NwNeVg2g+d7kMppVShtBpKKaVUoTRZKKWUKpQmC6WUUoXSZKGUUqpQmiyUUkoVSpOFUkqpQmmyUEopVShNFkoVAxH5i4iscW05+V6X6xHyqvTQQXlKFSPXWgg/G2PqOx2LUp7QkoVSxasFsM7pIJTylCYLpYpXS2C900Eo5SlNFkoVLy1ZqFJJk4VSxUtLFqpU0gZupYqJiAQAR4DqxpgMh8NRyiNaslCq+FyIXWlOE4UqdbRkoZRSqlBaslBKKVUoTRZKKaUKpclCKaVUoTRZKKWUKpQmC6WUUoXSZKGUUqpQmiyUUkoV6v8BZMToJ8m2Iz8AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "tmax = 10\n",
    "tau = 0.5           # transmission rate\n",
    "gamma = 1.0    # recovery rate\n",
    "rho = 0.05      # random fraction initially infected\n",
    "\n",
    "t, S, I, R = EoN.fast_SIR(BA, tau, gamma, rho=rho, tmax = tmax)\n",
    "plt.plot(t, S, color = 'r', label='S')\n",
    "plt.plot(t, I, color = 'b', label='I')\n",
    "plt.plot(t, R, color = 'g', label='R')\n",
    "plt.legend(loc=0)\n",
    "plt.xlabel(\"$T$\")\n",
    "plt.ylabel(\"individuals\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "caa992ff",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
